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ABSTRACT 


The  Structural  Synthesis  Transformation  is  used  to  conduct 
structural  system  identification  in  the  frequency  domain.  For 
spatially  complete  cases  where  each  of  the  frequency  response 
functions  at  every  degree  of  freedom  of  each  of  the  coordinates  of 
the  modeled  system  are  available  it  is  shown  that  the  theory  exactly 
identifies  all  modeling  errors.  For  spatially  incomplete  cases  where 
the  frequency  response  functions  are  available  only  at  a  proper 
subset  of  the  degrees  of  freedom  of  the  finite  element  model,  single 
mode  solutions  are  computed  over  intervals  about  the  modes  of  the 
experimental  system  using  matrices  and  complex  valued  line 
integrals.  Methods  of  forming  multiple  mode  solutions  from  the 
single  mode  solutions  are  explored. 
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I.  INTRODUCTION 


The  Finite  Element  (FE)  method  is  a  proven  tool  for  modeling  structural  dynamic 
systems.  As  the  complexity  of  the  system  increases  the  FE  model  may  not  accurately 
reflect  the  dynamic  behavior  of  the  system.  To  determine  the  extent  to  which  the  FE 
model  accurately  describes  the  physical  system,  a  comparison  of  the  dynamic  response  or 
modal  parameters  of  the  system  as  predicted  by  the  FE  model  and  the  re^onse  or  modal 
parameters  of  the  physical  system  as  determined  by  measurements  of  the  dynamic 
response  of  the  ^stem  can  be  made.  Such  a  comparison  can  easily  point  out  the 
differences  in  the  dynamic  behavior  of  the  FE  model  and  the  physical  system  but  fail  to 
provide  the  necessary  corrections  to  the  FE  model  that  will  provide  a  more  accurate  FE 
model  prediction  of  the  dynamic  behavior  of  the  physical  system. 

Structural  ^stem  identification  refers  to  procedures  used  to  identify  finit«:>  element 
modeling  errors  using  dynamic  test  data.  Localization  is  the  process  of  identifying  those 
degrees  of  freedom  (DOF)  of  the  FE  model  whose  in^edance  differs  from  that  of  the 
physical  system.  We  refer  to  this  set  of  DOF  as  the  error  set.  Identification  is  the  process 
of  finding  matrices  hK,  AM,  and  AC  that  are  corrections  to  the  FE  model  stiffiiess,  mass 
and  damping  matrices.  When  corrections  are  installed,  the  frequency  response  of  the 
corrected  FE  model  better  predicts  the  frequency  response  of  the  experimental  system, 
and  hence  the  modal  parameters  of  the  corrected  FE  model  better  predicts  those  of  the 
e}q)erimental  system. 

Structural  S5^em  identification  can  be  conducted  using  either  modal  or  frequency 
domain  methods.  This  thesis  investigates  a  frequency  domain  method  based  on  the 
structural  synthesis  transformation  (SST)  as  outlined  in  Reference  (1).  When  measured 
dynamic  response  data  is  available  fr)r  each  DOF  of  each  coordinate  of  the  FE  model  the 
identification  is  termed  spatially  complete  and  the  SST  can  be  used  to  exactly  determine 
which  elements  of  the  FE  model  are  in  error.  Additionally,  the  SST  can  be  used  to 
confute  correction  matrices  AM,  AK,  and  AC  that  can  be  used  to  correct  the  mass, 
str&ess,  and  damping  matrices  of  the  FE  moi^l.  The  dynamic  response  of  the  corrected 
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FE  model  exactly  matches  that  of  the  physical  system.  In  the  more  frequent  case  that 
dynamic  response  data  is  not  available  at  every  DOF  of  the  FE  model,  the  identification  is 
termed  spatially  incomplete.  In  this  case  the  SST  provides  a  frequency  dependent 
solution. 

SST-based  structural  system  identification  uses  the  frequency  reqjonse  function 
(FRF)  of  the  structure  under  consideratioiL  The  FRF  is  a  quantitative  description  of  the 
dynamic  behavior  of  the  structure.  To  obtain  the  FRF,  known  harmonic  excitation  forces 
are  applied  and  the  resulting  harmonic  response  of  the  structure  measured.  The  ratio  of 
excitation  forces  to  response  at  a  coordinate  evaluated  at  each  frequency  in  a  specified 
bandwidth  defines  the  FRF. 
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n.  THEORY 


A.  IMPEDANCE  DESCRIPTION 

The  impedance  model  of  a  given  physical  system  can  be  defined  by  the  relationship 
of  structural  displacement  with  respect  to  an  applied  force, 


(2.1a) 


The  harmonic  force  and  response  vectors  are  denoted  by  "f '  and  "x"  respectively. 
These  vectors  and  the  impedance  matrix,  Z,  are  in  general  conq)lex- valued  and  fi’equency 
dependent.  Subscripts  "i"  and  "c"  denote  non-error  and  error  DOF  respectively.  The 
superscr^t  "a"  denotes  quantities  calculated  firom  a  FE  (analytic)  model.  If  the  values 
were  obtained  firom  ejiperimental  test  data,  the  superscript  would  be  "x".  Thus,  for  the 
ejqjerimental  model  the  impedance  relationship  would  be  given  by: 


In  pratice,  the  elements  of  the  experimental  in^dance  matrix  of  Equation  (2.1b) 
are  unmeasurable  quanties.  To  see  this  we  e^qsand  row  j  of  Equation  (2.1b)  to  obtain 

fj  =  +  4^2+*  •  •+2>„  (2.2) 

To  measure  an  element  Zyj  of  Z“,  requires  that  we  impose  a  unit  displacement  at 

coordinate  xi  while  physically  holding  all  other  other  coordinates  at  zero  displacements. 
This  is  not  physically  possible.  Assuming,  for  purpose  of  definition,  the  availability  of  the 
experimental  impedance  matrix,  the  quantitative  difference  between  the  analytical  and 
experimental  systems,  as  a  function  of  frequency,  is  described  by  the  impedance  error 
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matrix.  It  is  defined  by  the  difference  between  the  analytical  and  experimental  impedance 
matrices.  The  error  impedance  matrix  relationship  is  defined  as: 


’0  0  ■ 

1 — 

0  AZ(Q)_ 

ya  ya 

\_^ci  ^ccj 

yx  yx 
^cc 


(2.3) 


B.  STRUCTURAL  SYNTHESIS  TRANSFORMATION 

Since  the  ejqperimental  in^edance  matrix,  2*,  is  in  general  unavailable,  firequency 
domain  structural  synthesis  is  used  to  identify  the  impedance  error  matrix  nging  FRF  data 
exclusively.  A  structural  synthesis  transformation  is  constructed  fi-om  AZ  of  Equation  (2.3) 
which  encompasses  the  FE  model  errors.  This  transformation  is  applied  to  the  finite 
element  model  to  produce  an  e3q)erimental  system  FRF. 

The  FRF  relates  structural  response  to  applied  excitation.  Given  a  FE  model  with 
inqjedance  matrix,  2" ,  the  FRF,  If,  is  the  matrix  inverse  of  the  impedance  matrix  2“  of 
equation  (2.1a).  For  a  static  system  (Q=0)  the  FRF  is  simply  the  flexibility  matrix  (inverse 
stififiiess  matrix).  Using  the  notation  of  Equation  (2.1)  we  may  partition  If  as  follows: 


(2.4) 


Generalfy,  the  "c"  response  coordinates  e3q)erience  applied  forces  due  to  both  error 
inpedances  and  externally  applied  forces,  whereas  "i"  response  coordinates  experience 
only  externally  applied  forces,  such  that, 

fc=fT^f^  (2.5a) 


and 

f>=fr 


(2.5b) 


.  Expanding  Equation  (2.4)  and  substituting  the  relations  of  Equation  (2.5)  yeilds 
the  following  relationships: 

X,  =  HtfT  +  K/r  +  Kf^  (2.6a) 

X,  =  Hlfr  +  KfT  +  (2.6b) 


If  we  include  a  copy  of  Equation  (2.6b),  in  expanded  matrix  notation.  Equation 
(2.6)  reflects  the  three  harmonic  excitation  terms  to  be  considered,  i.e.. 


K 

K 

fT' 

►  = 

K 

K  ■ 

fT- 

K 

K 

rbZ 

Jc 

V  J 

Reqwnse  coordinates  "c"  and  "i",  which  are  due  to  external  forces,  will  hereafter 
be  referred  to  as  "e"  coordinates,  denoting  their  dependence  on  external  force  excitation. 
Consequently,  the  three  excitation  forces  are  condensed  into  two  under  the  identity: 


{/.}=[[/“]  [r]f 


(2.8a) 


U}  =  [ff] 


and  Equation  (2.7)  reduces  to 

'K  iClfXl 


(2.8b) 


(2.9) 


Equation  (2.3)  shows  that  the  impedance  error  is  defined  as  the  difference  between 
the  analytic  and  ejqjerimental  inqiedance  models.  Hence,  a  transformation  is  required 
which  uses  the  FRF  relationship  of  Equation  (2.9)  to  generate  a  similar  relationship  for  the 


5 


e3q)erimental  system.  The  impedance  error  AZ  provides  the  basis  by  which  this 
transformation  is  developed. 


The  inq)edance  error  matrix,  AZ(Q),  is  a  fimction  of  frequency  and  satisfies: 


{/4  =  -[Mn)]k} 


(2.10) 


where 


[az(q)] = [[A^:]  +  yQ[AC]  -  q^[am]] 


(2.11) 


for  Q  the  forcing  frequency,  j  =  and  AAT,  AC,  and  AM  stiflfiiess,  dancing  and  mass 
error  matrices  con^arable  to  those  of  the  finite  element  formulation.  The  minus  sign  in 
Equation  (2.10)  reflects  that  the  reaction  forces  imposed  by  inqjedance  errors  on  the 
baseline  model  are  being  considered.  Substituting  the  relationship 

{3=[o  -3K1 


into  Equation  (2.9)  yields: 


K  KU  0  irx 

HI  0  -AZU. 


(2.13) 


simplifying  we  get: 


HI  -^:Azir/, 


H:,  -/CAZjk 


(2.14) 


Expanding  Equation  (2.14)  into  two  equations  and  using  to  denote  a 
synthesized  modified  response  results  in 


(2.15a) 


x:=Hlf,-HlhZx- 


eeJe 


(2.15b) 
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Rearranging  Equation  (2. 15b)  produces 

\l  +  =  HU. 

(2.16a) 

Using  the  property  of  the  frequency  response  function, 

=  HU 

(2.16b) 

[I+H°^mYc=x^ 

(2.16c) 

x:  =[/  +  i/“AZ]‘*jc, 

(2.16d) 

Introducing  Equation  (2.16)  into  Equation  (2.15b)  results  in 

x;  =  HUe  -  ^CAZ[/ + Ki^Yx, 

(2.17a) 

= HU  -  ka2{i+h:,azYhu 

(2.17b) 

Once  again  we  recall  the  property  of  a  FRF, 

= HU. 

(2.18a) 

Combining  Equations  (2.17b)  and  (2.18a)  yields 

-H°,A2\i+H^AzYh°, 

(2.18b) 

Noting  that 

[/+  fCAZ)"'  =  [(aZ-‘  +  /C)A2f ' 

(2.19a) 

and  applying  the  matrix  property 

(2.19b) 

we  get  that 

a;  =  «;-«■![  AZ-'  +  fljf/i; 

(2.20a) 
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Replacing  the  superscript  which  denotes  the  structures's  synthesized  coupled 
response,  with  the  superscript  "x"  to  indicate  the  test  ^stem  response  we  arrive  at 

HI- h:\6Z-' (2.20b) 

In  full  matrix  notation  we  have 


HI 

Hli 

H^cc 

rrd  rja 
Ilci  ^cc 


K 

.K 


[az-' + 


HI 

m 


(2.20c) 


Equation  (2.20b)  is  the  structtiral  ^thesis  transformation  equation.  When  the 
experimental  system  FRF,  If  ,  is  available  the  SST  can  be  used  to  identify  a  frequency 
dependent  impedance  error  matrix  [AZ(Q)].  Additionally,  using  Equation  (2.9),  [AZ(Q)] 
can  be  decomposed  into  constitutent  stiffiiess,  mass,  and  damping  errors. 


C.  FREQUENCY  DOMAIN  LOCALIZATION 
We  may  rewrite  Equation  (2.20b)  as 

AH^  =  (2.21) 

where 

^ee=Hte-H^ee  (2.22a) 

and 

£>  =  [AZ'‘+i/“]  (2.22b) 

We  define  the  localization  matrix  L  as 

L  =  (2.23) 

using  the  expression  of  Equation  (2.22a)  in  Equation  (2.23)  we  can  rewrite  L  as 

Z  =  (2.24a) 

Expanding  the  "e"  coordinate  set  into  error  and  non  error  coordinates  we  get 
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(2.24b) 


L  = 


Z 


ya 

ya 


7" 

^ic 


Noting  the  the  frequency  response  matrix  is  the  inverse  of  the  impedance  matrix,  i.e., 


1 

'K 

K 

■/ 

O' 

[z: 

I _ 

0 

(2.25a) 


all  mixed  product  coordinates  of  Equation  (2.24b)  must  be  zero  and  L  simplifies  to 


'0  0  ■ 
0  £>■' 


@  Q=Qi 


(2.25b) 


D.  ERROR  IMPEDANCE 

We  can  solve  Equation  (2.20b)  for  the  in5)edance  error  [AZ].  From  Equation 
(2.25b)  terms  of  Equation  (2.20b)  associated  with  non  error  coordinates  may  be  assumed 
to  be  zero.  We  get  the  ft)llowing  fijrm  of  the  impedance  error  matrix 

[AZl  =  ([ffi]-[/fj])''  @n=n,  (2,26a) 

where 

=  (2.26b) 

Let  S={ni,  Q2,  ... ,  fin}  be  a  set  of  frequencies  where  Qi  <  n2<...<Qn-i<Qn.  If 
for  each  i=  2,3, ...,  n-1,  we  apply  Equation  (2.10)  at  each  of  the  frequencies  Qm,  Qi ,  and 
Qj+i  and  assemble  the  resulting  three  equations  into  a  system  of  three  equations  in  three 
unknowns,  we  get  the  matrix  equation 


J^i-/ 

= 

I  -qU 

jClJ 

AM, 

I  -Ql,I 

Equation  (2.27)  can  be  used  to  decompose  the  frequency  dependent  in:^>edance 
error  into  constitutent  stiffiiess,  mass,  and  dancing  error  matrices,  AK,  AM,  and  AC  at  the 
frequencies  Qj,  i— 2,3,...,n-l.  These  constitutents  matrices  are  in  general  frequency 
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AK^(Q,)- 

dependent.  Denoting  the  solution  of  Equation  (2.27)  by  ]  AA<,(Q;)  ►  for  i=2,3,...,n-l  we 

AC(n,) 

obtain  for  each  Qi  i=2,3,...4i-l,  error  stiffoess,  mass,  and  damping  matrices  AKc(Qi), 
AMc(f2i),  and  ACc(f2i)  such  that 

Az,(n,) = Ax,(n,)  -  Q?Aiw;(Q,) +mac,(q,)  (2.28) 

The  matrices  AKc(Oi),  AM.c(Q^,  and  ACcfjQi)  can  be  used  to  numerically  correct 
the  stififiiess,  mass,  and  damping  matrices  of  the  FE  model  at  the  frequencies  Qi 
i~2,3,...,n-l  in  that  the  ITRF  of  the  corrected  FE  model  at  the  frequencies  Q,  approximates 
the  experimental  system  FRF  at  Qi.  To  express  this  symboKcaUy,  ifK°,  Kf,  and  C  are  the 
stiffiiess,  mass,  and  damping  matrices  of  the  FE  model  and  If  is  the  FRF  matrix  of  the 
experimental  system  at  Qi 

(Q,)  « [r  +  A^:,(Q,)  -  Q2  +  AM,(Q,)) + yQ  .(c°  +  AC,(Q,)))''  (2.29) 

where  the  "c"  subscripted  matrices  are  added  at  the  corresponding  error  set  coordinates  of 
the  "a"  superscripted  matrices. 

For  a  general  set  of  frequencies,  S={Qj,  Qj, ... ,  Q„},  we  can  form  the  system  of  n 
equations  in  three  unknowns  given  by: 

'AZ,(Qj)1  P  -Q2/  jQ^f 

i  :  :  :  Ta/s:/ 

AZ,(Q,)  =  /  -Qf/  yQj  (2.30) 

:  i  i  : 

AZ,(Q„)J  [l  -all  JQJ 

rA/rp  rA^,(S)' 

The  solution  AM^  —  AM^(a)  of  Equation  (2.30)  represents  error  stiffiiess, 

_ACJ  LaQ(5)_ 

mass,  and  danqiing  matrices  which  best  corrects  the  FE  model  in  a  least  squares  sense. 
Equations  (2.27)  and  (2.30)  are  fundamental  to  all  that  follows. 
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m.  SPATIALLY  COMPLETE  STRUCTURAL  IDENTIFICATION 


To  illustrate  the  principles  of  SST  based  frequency  domain  structural  identification 
we  will  use  the  FE  model  of  a  free-free  beam.  To  simulate  the  experimental  system  we 
impose  a  25%  addition  to  the  mass  and  stiffiiess  of  elements  3  and  4  of  a  10  element  FE 
model.  Figure  3-1  shows  the  finite  element  model  of  the  beam  and  the  spatially  complete 
experimental  system  that  results  from  imposing  the  mass  and  stiflfiiess  additions  at  element 
3  and  4  of  the  FE  model.  Table  3-1  shows  the  system  frequencies  of  the  analytic  and 
experimental  systems. 


EXPERIMENTAL  SYSTEM 


MASS  OF  ELEMENT  3  =  1.954X10^  LBM 
MASS  OF  ELEMENT  4  » 1.954X10*-6  LBM 
STIFFNESS  OF  ELEMENT  3  =  417.88  LBF/IN 
STIFFNESS  OF  ELEMENT  4  »  417.88  LBF/IN 


Figure  3- 1  Spatially  complete  Analytic  and  Experimental  Systems. 
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We  will  denote  by  M” ,  K‘‘ ,  and  €“  the  mass  stiffiiess  and  damping  matrices  of 
the  FE  model  and  by  ,  and  the  mass  stifi&iess  and  damping  matrices  of  the 

experimental  system.  The  impedance  matrices  of  the  analytic  and  simulated  experimental 
systems  are  given  by; 
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Z‘’(Q)  =  K‘’+  jac  - 


(3.1a) 


Z=^{a)  =  K^  +jQ.C 


(3.1b) 


The  FRF  matrices  of  the  analytic  and  simulated  e}q)erimental  systems  are  given  by: 


/f“(Q)  =  Z“(Q)"^ 

(3.2a) 

/f"(fi)  =  Z^(Q)'’ 

(3.2b) 

Figure  3-2  shows  a  comparison  of  the  driving  point  FRF  at  DOF  1  of  the  analytic 
and  experimental  system. 
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Figure  3-  2  Analytic  FRF  vs  simulated  Experimental  FRF. 
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For  a  given  frequency  Qo,  using  Equation  (2.24a),  we  can  form  the  localization 
matrix  L.  Plotting  the  diagonal  elements  of  L  versus  the  associated  DOF  we  obtain  the 
plot  shown  in  Figure  3.3.  For  each  frequency,  Q,  in  a  given  frequency  range,  we  can 
con:q)ute  the  diagonal  of  L  at  Q  .  Assembling  the  diagonals  over  the  frequency  range  of 
our  system  into  a  rectangular  matrix  and  performing  a  MATLAB  mesh  plot  of  the 
resulting  rectangular  matrix  we  obtain  the  surface  plot  shown  in  figure  3.4.  Figure  3.4 
shows  the  frequency  dependency  of  the  localization  matrix  diagonals.  As  our  system  is 
spatially  conplete  Equation  (2.25)  forces  all  non  error  coordinates  to  be  zero.  From 
Figure  3-3  we  can  determine  that  the  locations  of  the  nonzero  diagonal  values  are  DOF  9, 
10,  11,  and  12.  We  denote  by  the  set  of  DOF  for  which  L(i,i)  is  non  zero.  For  our 

beam  system  =  {9,10,11,12} .  Figure  3.5  shows  the  frequency  dependence  of  typical 
error  and  non  error  set  diagonal  elements  of  the  localization  matrix  L. 

Using  the  set,  Cerr»  which  results  from  the  localb^tion,  we  can  perform  a 
partitioning  of  the  FRF  matrix  as  described  by  Equation(2.4).  We  can  now  apply 
Equation  (2.26)  to  compute  AZ  as  a  matrix  fimction  of  frequency  over  the  frequency 
range  of  our  system.  We  then  use  Equation  (2.27)  to  decompose  A7,  into  its  constituent 
components  AK,  AM,  and  AC.  We  get  exact  solutions  of  error  stiflSiess,  mass,  and 
damping  as  shown  in  figures  3-7,  3-8  and  3-9.  The  MATLAB  Routine  SST.M  of 
Appendbc  A  can  be  used  to  accomplish  the  above  steps. 

For  each  Q  in  the  frequency  range  of  our  system  we  can  form  the  sum: 

=  Z“  (Q)  +  AZ(Q)  (3.3) 

We  refer  to 

^corr  ~  (3  >4) 

as  the  corrected  FRF  of  the  analytic  model.  Hcorr  is  the  FRF  of  the  model  that  results  from 
installing  the  corrections  as  identified  by  the  Equation  (2.26a).  Figure  3-10  shows  a 
comparison  plot  of  the  FRF  of  the  corrected  model  and  the  FRF  of  the  experimental 
system.  Figure  3-10  clearly  shows  the  exactness  of  the  SST  solution  in  the  case  of  a 
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spatialfy  complete  system.  The  experimental  and  corrected  model  FRFs  are  identical  to 
within  plot  resolution. 

At  each  frequency,  Q,  in  the  frequency  range  of  our  system  we  can  form  the  sum 
(Q)  =  (Q)  +  A^(Q)  +  jOACifi)  -  (3.5) 

We  refer  to: 

(3.6) 

as  the  constituent  corrected  FRF.  Figure  3-11  is  a  conparison  plot  of  the  constitutent 
corrected  FRF  and  experimental  for  our  system  The  offeet  between  the  two  plots  is 
exactfy  equal  to  the  sanpling  frequency  AQ. 
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Figure  3-  3  Spatially  complete  localization  matrix  diagonal  at  f2=352  Hz. 
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Non-Error  Coord  L(14,14)  Error  Coord  L(9,9) 


Figure  3-  5  Frequency  dependence  of  spatially  complete  localization  matrix  error  set 
DOF.  Units  are  Ibf'in  (LoglO  o:^- 
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Z(7,7)  lbf/in(LoglOof) 


Figure  3-  6  Spatially  complete  Analytic  Impedance  vs  Experimental  impedance. 
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Figure  3-  7  Computed  stiflfiiess  vs  true  stiffiiess  for  spatially  complete  beam.  Plots  are 
identical  within  plot  resolution. 
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Figure  3-  8  Computed  mass  vs  true  mass  for  spatially  complete  beam.  Plots  are  identical 
within  plot  resolution. 
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Figure  3-  9  Computed  damping  vs  true  damping  for  spatially  complete  beam.  Plots  are 
identical  within  plot  resolution. 
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Figure  3-  10  Spatially  conplete  Ejqjerimental  FRF  vs  AZ  corrected  FRF.  Plots  are 
identical  within  plot  resolution. 
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Figure  3-11  Spatially  complete  e)q)erimental  FRF  vs  AK,  AM,  and  AC  corrected  FRF. 
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IV.  SPATIALLY  INCOMPLETE  STRUCTURAL 
IDENTIFICATION 

A.  GENERAL  DESCRIPTION 

To  illustrate  SST  based  frequency  domain  structural  identification  when  the 
physical  system  under  consideration  is  spatially  incomplete  we  will  again  use  the  FE  model 
of  a  free-free  beam.  We  simulate  the  ejqjerimental  system  by  imposing  a  25%  addition  to 
the  mass  and  stiffiiess  of  elements  3  and  4  of  a  10  element  FE  model.  Figure  4-1  shows  the 
FE  modeled  beam  and  the  spatially  incomplete  system  that  results  if  FRF  data  is  available 
only  at  the  displacement  DOF  of  the  FE  system. 
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Figure  4- 1  Analytic  and  experimental  spatially  incomplete  systems. 
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To  obtain  a  simulated  FRF  for  our  spatially  incomplete  beam  we  denote  by 
M° ,  K° ,  and  C"  the  mass  stififiiess  and  damping  matrices  of  the  FE  model  and  by  M' , 
,  and  the  mass  stifl&iess  and  damping  matrices  of  the  experimental  ^stem.  Af' , 
AT' ,  and  C'  are  obtained  by  imposing  15%  mass  and  stifl&iess  and  a  15%  mass  additions 
to  the  elemental  matrices  of  element  5  and  6  respectively  of  the  FE  model.  The  impedance 
matrix  of  the  simulated  spatially  complete  beam  is  given  by 

Z"(Q)  =  jClC  -  (4-1) 

The  FRF  matrix  of  the  simulated  spatially  conq)lete  system  is  given  by 

H\a)  =  .  (4-2) 

We  introduce  the  terminology  Analysis  set  and  Omitted  set  where  the  Analysis  set 
(A-set)  is  that  set  of  DOF  for  which  experimental  FRF  data  is  available  and  the  Omitted 
set  (0-set)  is  that  set  of  DOF  of  the  ej^erimental  system  for  which  experimental  FRF  data 
is  imavailable.  For  our  simulated  ejqjerimental  system  the  A-set  consist  of  the  odd 
numbered  translational  DOF  and  the  0-set  consists  of  the  even  numbered  rotational  DOF, 
i.e., 

=  {1,3,5, ...,21}  (4.3a) 

(2,4,6,..., 22}  (4.3b) 

We  obtain  the  simulated  FRF  of  our  spatially  incomplete  beam  by  physically 
extracting  the  rows  and  columns  of  the  simulated  spatially  complete  matrix,  If,  for  which 
FRF  data  would  be  available.  In  our  system  these  are  the  rotational  DOF  and  all  even 
numbered  rows  and  columns  are  omitted  from  the  simulated  spatially  cor^plete  FRF  to 

obtain  a  spatialfy  inconqjlete  FRF  which  we  will  denote  hyW . 

For  a  fixed  Qo,  is  a  square  matrix  of  size  (l&igth(A-set))  by  (length(^- 

sef)).  For  our  FE  model  as  currently  defined.  If,  is  of  size  (number  of  DOF)  by  (number 

of  DOF)  and  is,  as  is  true  of  most  real  world  cases,  of  larger  size  than  W .  In  order  to 
employ  the  structural  synthesis  transformation  we  need  to  reduce  the  size  of  If  to  that  of 
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H  .  To  this  end  we  will  consider  two  reduction  methods,  FRF  matrix  extraction  pief  1], 
and  the  Improved  Reduced  System  as  given  in  [Ref.  2]. 


B.  EXTRACTION  REDUCTION  METHOD 

In  order  to  reduce  If  by  the  extraction  method  we  simply  extract  from  the  full 
order  H  those  rows  and  column  which  correspond  to  A-set  coordinates.  Partitioning  the 
impedance  and  fuU  FRF  matrices  of  our  analytical  system 


(4.4a) 


(4.4b) 


(4.5) 


(4.6) 


C.  O-SET  SYSTEM 

Equation  (4.6)  relates  the  extracted  reduced  order  FRF  of  the  analytic  model  to  the 
impedance  of  the  full  order  order  model  Taking  advantage  of  the  identity 

[^ooY  =  {dt\z^]f'adj[z^]  (4.7) 

and  replacing  Zoo  by  Koo+j£2Coo-jff Moo  (where  det(*)  and  adj(*)  represent  the 
determinant  and  adjount  respectively),  we  see  that  an  element  is  large  for  those 

frequencies  Qo  where  Qo  is  an  eigenvalue  of  the  0-set  system,  the  0-set  system  being  that 
FE  model  having  stififiiess,  mass  and  danqiing  matrices  Koo,  Moo,  and  Coo  respectively. 
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D.  IMPROVED  REDUCTION  SYSTEM 


To  xise  the  in:5)roved  reduction  method  (IRS)  we  first  partition  2“  using  the  A  and 
O  sets,  then  adjust  Equation  (2.1)  to  reflect  this  new  coordinate  system  obtaining 


II 

^aa  ^ao 

(4.8) 

E^qjanding  Equation  (4.8)  into  two  equations  yields 

fa  =  +  KoXo 

(4.9a) 

fo- 

(4.9b) 

0-set  coordinates  are  not  associated  with  FRF  data  measurement  locations  on  the 
physical  structure,  therefore  the  forcing  fimction  at  O-set  coordinates  can  be  set  to  zero. 
Making  these  substitution  into  Equation  (4.9b)  and  solving  for  the  generalized  structural 
response  coordinates  leads  to: 


=  —7-1 7"  V 


Substituting  these  results  into  Equation  (4.8)  yields, 


f/.l 

l- 

^aa 

Zi- 

/  ■ 

loj 

r 

_oa 

1 

N 

^oo 

hence 


(4.10a) 

(4.10b) 

(4.11) 

(4.12) 


When  Q=0  the  Equation  (4.10a)  yields  the  static  reduction  relationship  between 
omitted  and  retained  coordinates  and  is  given  by 


(4.13) 
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The  IRS  relationship  is  given  by 

M  =  \-K-^K^  +  TM-Z,K„\{x.}  (4.14) 

where 

T  =  K2M^-K-^MJC2K„  (4.15) 

and  Kstat  and  Msm  are  the  statically  reduced  [Ref.  2, 3]  stifiSiess  and  mass  matrices. 

Unlike  the  spatially  complete  case  where  we  only  had  to  consider  two  system  (the 
analytic  and  the  experimental  systems),  in  the  case  of  a  spatially  incomplete  system  there 
are  actually  five  systems  with  which  we  must  concern  ourselves  [Ref  4]:  the  anafytic 
system,  the  experimental  system,  the  reduced  analytic  system  that  results  from  conducting 
dynamic  reduction  on  the  mass  and  stiffiiess  matrices  of  the  analytic  system,  anH  the 
omitted  systems  of  both  the  analytic  and  the  experimental  system  Table  4-1  shows  the 
frequencies  of  each  the  first  four  of  these  systems. 

Figure  4-2  shows  a  comparision  of  the  analytic  and  experimental  FRF  of  our 
spatially  incomplete  beam  We  see  that  by  using  IRS  reduction  of  the  analytic  system  those 
modes  above  approximately  1000  Hz  i.e.,  those  modes  associated  with  the  reduced  out 
rotations,  are  not  present.  Figure  4-3  shows  a  similar  comparision  where  we  have  used 
extraction  reduction  and  the  higher  modes  are  present.  In  all  that  follows  we  shall  use  IRS 
reduction  of  our  analytical  systems.  Figure  4-4  shows  the  localization  matrix  diagonal  at 
^“196.1  Hz  while  Figure  4-5  shows  the  frequency  dependency  of  the  localization 
diagonals  over  the  frequency  range  our  our  spatially  incomplete  beam  For  a  spatially 
incomplete  system  the  determination  of  the  locations  of  the  error  coordinates  is  not  a 
clearly  defined  task.  We  shall  not  discuss  the  problem  of  actually  determining  the  exact 
error  set  in  the  spatially  inconq)lete  case  and  will  use  our  knowledge  of  the  true  location  of 
the  error  coordinates  to  aid  in  our  localization.  As  we  shall  make  use  of  the  concept  of  the 
size  of  the  error  set  again,  we  will  sinqjly  note  that  using  a  reduction  method  like  IRS 
causes  errors  to  be  'smeared'  in  the  reduced  anatytic  model.  Table  4-2  shows  the  A-set  and 
0-set  of  our  simulated  beam  along  with  the  locations  of  the  true  errors  as  well  as  the 
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computed  C-set  at  Q=196  Hz.  The  computed  C-set  is  the  set  of  all  DOF  having  a  diagonal 
entry  whose  absolute  value  exceeds  a  given  tolerance. 


Mode 

Anal  (Hz) 

Exp  (Hz) 

Reduced 

(Hz) 

0  set  (Hz) 

1 

25.41 

25.83 

25.41 

1244.28 

2 

70.07 

70.85 

70.07 

1287.94 

3 

137.45 

137.93 

137.45 

1417.14 

4 

227.55 

227.37 

227.55 

1629.59 

5 

340.84 

341.14 

340.84 

1928.6 

6 

478.0 

477.97 

478.0 

2327.83 

7 

639.9 

639.92 

6403 

2853.66 

8 

826.5 

826.24 

829.7 

3541 

9 

1026.8 

1025.99 

1029.74 

4399.32 

10 

1363.62 

1364.51 

5280.86 

11 

1640.62 

1640.71 

5701.99 

12 

1981.19 

1984.61 

13 

2381.4 

2380.32 

14 

2850.01 

2847.38 

15 

3397.39 

3397.74 

26 

4027.82 

4026.28 

17 

4720.15 

4738.64 

18 

5376.83 

5380.09 

19 

6794.91 

6787.09 

20 

6807.52 

6804.12 

Table  4-1  Analytic,  Experimental,  Reduced,  and  Omitted  System  Frequencies  of  Spatially 
Inconq)lete  Beam. 
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A  Set 

1 

3 

5 

7 

9 

11 

13 

15 

17 

19 

21 

OSet 

2 

H 

6 

8 

10 

12 

14 

16 

18 

20 

22 

True  C  Set 

5 

6 

B 

8 

9 

10 

B 

B 

B 

B 

Computed  C 

Set 

1 

3 

5 

7 

9 

11 

13 

1 

1 

1 

1 

Table  4-2  Analytic  set,  Omitted  set.  Computed  C  set,  and  True  C  set  DOF  for  spatially 
incomplete  beam. 
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Figure  4-  2  Analytic  FRF  vs  ejqserimental  FRF  for  spatially  incomplete  beam  using  IRS 
reduction. 


Figure  4-  3  Analytic  FRF  vs  experimental  FRF  for  spatially  incomplete  beam  iistng 
extraction  reduction. 
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L(DOF)  Ibf/in 


ASET  DOF 


Figure  4-  4  Localization  matrix  diagonal  at  0=196.1  Hz  for  spatially  incomplete  beam. 
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Figure  4-  5  Frequency  dependence  of  localization  matrix  diagonals  for  spatially 
incomplete  beam. 
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Non-Error  Coord  L(15,15)  Error  Coord  L(1 1,1 1) 
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Figure  4-  6  Frequency  dependence  of  error  and  non  error  set  localization  matrix 
diagonals  for  spatially  incon:5)lete  beam.  Units  are  Min  (LoglO  o^. 
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Omega''2  (Hz) 


Figure  4-  7  Analytic  vs  experimental  impedance  for  spatially  incomplete  beam. 


38 


DM(7,7)  lbf-sec'^2/m  (LoglO  of) 


Figure  4-  9  Computed  Mass  vs  True  Mass  for  spatially  incomplete  beam. 
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DC(7,7)  Ibf-sec/in  (LoglO  of) 


Omega  (Hz) 


Figure  4- 10  Computed  Damping  vs  True  Damping  for  spatially  incomplete  beam. 
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H(7,7)  iii/lbf(LoglOof) 


Figure  4- 11  AZ  Corrected  FRF  vs  e3q)erimental  FRF  for  spatially  incomplete  beam.  Plots 
are  identical  to  within  plot  resolution. 


Figure  4-  12  AK/AM/AC  Corrected  FRF  vs  experimental  FRF  for  spatially  incortplete 
beam.  The  oflfeet  of  the  two  plots  is  equal  to  the  sampling  frequency  AQ. 
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V.  SINGLE  MODE  SOLUTIONS 


A.  SINGLE  MODE  MATRIX  SOLUTIONS 

In  the  case  of  a  spatially  conq)lete  system  we  have  seen  that  the  SST  jdelded 
frequency  independent  error  matrices  AK,  AM,  and  AC  that  could  be  used  to  correct  the 
stiffiiess,  mass,  and  damping  matices  of  the  FE  model  in  such  a  manner  that  the  FRF  of  the 
corrected  FE  model  was  exactly  equal  to  that  of  the  e}q)erimental  systenx  In  the  case  of  a 
spatially  incomplete  system  the  constituent  solutions  AK,  AM,  and  AC  given  in  chapters 
III  and  IV  were  in  general  frequency  dependent  solutions  which  serve  only  as  corrections 
to  the  reduced  FE  model.  Ideally  we  seek  frequency  independent  solutions  which  are 
corrections  to  the  full  FE  model.  For  now  we  shall  only  deal  with  the  simpler  problem  of 
trying  to  find  frequency  independent  solutions  which  serve  as  corrections  to  the  reduced 
FE  model. 


We  shall  employ  Equation  (2.30).  For  a  given  experimental  system  mode,  cOn, 
consider  a  frequency  bandwidth  [Qi,  QJ  such  that  (»n€[Qi,  QJ-  Let  S={Qi  Q2,...,Om}  be 
a  frequency  sampling  of  the  bandwidth  [Qi,  QJ,  i-e.,  Qi<Qi^u  for  each  QigS.  For  each 
QigS  we  can  form  the  error  impedance  matrix  AZ(Q).  We  apply  Equation  (2.30)  to  the 
partition  yeilding  the  following  system  of  m  equations  in  three  imknowns: 


■AZ.(n,)' 

7  jQJ 

*  .  .  r 

Az.(a,) 

I  -Q.]I  jQJ 

AZ,(n„) 

I  -Qil  jQJ 

ak;- 

am;- 

Ac;» 


(5.1) 


We  will  demonstrate  by  example  that  for  properly  chosen  bandwidths  [Qi,  QJ,  the 
solution  AK®",  AM^" ,  AC®"  of  Equation  (5.1)  approximately  corrects  the  Analytic  FRF 
over  the  frequency  bandwidth  [Qi,Qu],  i.e.,  if  H‘(0)  is  the  value  of  the  FRF  matrix  of  the 
experimental  system  at  a  frequency  Q,  //“(Q)the  value  of  the  FE  reduced  analytic  FRF 


44 


matrix  at  Cl,  and  H°^{Cl)  the  value  of  the  single  mode  corrected  reduced  Analytic  FRF 
matrix  at  Q  for  an  experimental  system  mode  tOn  where  Q  e  [Qijfiu]  then 


(5.2) 


We  will  refer  to  the  solution  AK^" ,  AM®" ,  AC®"  as  a  single  mode  solution  at  cOn. 

In  the  case  of  a  spatially  complete  system,  such  as  the  spatially  complete  beam  discussed  in 
chapter  III,  all  single  mode  solutions  are  foimd  to  be  identical  to  the  unique  frequency 
independent  solution  that  was  obtained  in  chapter  HI,  hence  the  corrections  are  valid 
throughout  the  frequency  range  of  the  experimental  system  and  the  left  hand  side  of 
Equation  5.2  is  zero.  For  a  spatially  incomplete  system  this  is  not  the  case. 

In  support  of  Equation  (5.2)  we  chose  a  1  Hz  frequency  bandwidth  centered  on 
mode  1  ((01=25.21  Hz)  of  our  spatially  incomplete  beam  as  defined  in  Chapter  III.  For  our 
frequency  sampling,  H,  we  will  use  a  sanqjling  frequency  of  .5  Hz  which  results  in  3  points 
which  are  equally  spaced  over  the  interval,  S={24.7178  Hz,  25.2178  Hz,  25.7178  Hz}. 
Solving  Equation  (5.1)  we  get 


ak:'  = 


am:'  = 


■  -3.91 

-13.60 

21.70  ■ 

-13.60 

170.5 

-187.5 

(5.3a) 

-21.70 

-187.5 

194.6  _ 

-0.0001 

0.0006 

-0.0005‘ 

0.0006 

-0.0033 

0.0026 

(5.3b) 

-0.0005 

-0.0026 

-0.0016 

Ac;- 


0  +  0.0270y  0- 0.229  ly  0 +  02277 J 
0-0.229iy  0  +  1.5078y  0-1.4627y 
0  +  02177 j  0-1.4627y  0  +  1.378iy 


(5.3c) 


Figure  5-1  shows  a  comparison  of  experimental,  uncorrected  and  mode  1 
corrected  FRFs  for  our  spatially  incomplete  beam  using  a  1  Hz  frequency  bandwidth 
centered  on  Mode  1  and  a  frequency  sampling  consisting  of  3  frequencies  equally  spaced 
over  the  bandwidth.  Figure  5-2  shows  the  results  of  including  mode  1  and  mode  2  in  the 
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frequency  bandwidth  [Qi,  QJ.  As  Figure  5-2  shows  Equation  5.1  is  not  as  accurate  when 
multiple  frequencies  are  included  in  the  bandwidth.  Figure  5-3  is  a  comparison  plot  of  the 
e3q)erimental  and  analytic  FRF  versus  the  single  mode  matrix  solution  corrected  FRFs  over 
25,  10,  and  1  Hz  bandwidths  using  a  3  point  frequency  sampling.  Figure  5-3  shows  that 
the  precedure  is  sensitive  to  the  size  of  the  bandwidth  over  which  the  solution  is 
computed,  i.e.,  better  accuracy  is  achieved  with  smaller  bandwidths.  Figure  5-4  is  a 
comparison  plot  of  the  experimental  and  uncorrected  analytic  FRFs  versus  the  single  mode 
matrix  solution  corrected  FRFs  computed  over  a  1  Hz  bandwidth  using  frequency 
samplings  of  200,  50,  10,  and  3  points  equally  spaced  over  the  bandwidth.  Figure  5-4 
shows  that  the  procedure  is  feirly  insensitive  to  sampling  frequency. 

Figure  5-5  is  a  comparison  plot  of  the  ejqjerimental  and  analytic  FRFs  versus  the 
single  mode  solution  corrected  FRF  for  a  1  Hz  bandwidth  using  a  3  point  frequency 
sampling  in  the  case  of  a  spatially  complete  beam.  As  Figure  5-5  shows  the  procedure  is 
exact  for  spatially  complete  systems. 
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H(7,7)  in/lbfaoglOof) 


H(7,7)  in/lbfaoglOof) 


Figure  5-  2  Experimental  and  imcorrected  analytic  FRFs  vs  single  mode  corrected  FRF 
using  a  15  point  frequency  sampling  of  a  bandwidth  that  includes  modes  1  and  2. 
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H(7,7)  iii/lbf(loglOof) 


Figure  5-  3  Spatially  incomplete  experimental  FRF  vs  single  mode  matrix  solutions  at 
mode  1  using  3  point  frequency  samplings  of  25,  10,  and  1  Hz  bandwidths. 
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Figure  5-  4  Spatially  incomplete  experimental  FRF  vs  single  mode  matrix  solutions  at 
mode  1  using  3, 10,  50,  and  200  point  frequency  samplings  of  a  1  Hz  bandwidth. 
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H(7,7)  in/lbf(loglOof) 


Figure  5-  5  Spatially  complete  eiqperimental  FRF  vs  single  mode  matrix  solutions  at 
mode  1  using  a  3  point  frequency  samplings  of  a  1  Hz  bandwidth. 
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B.  SINGLE  MODE  INTEGRAL  SOLUTIONS 

In  what  follows  we  wish  to  employ  the  integral  formulas  of  reference  (5)  to  obtain 
a  results  similar  to  that  of  Equation  (5.2).  To  accomplish  this  we  shall  need  to  recast  the 
impedance  equations  of  Chapter  n  in  terms  of  velocity.  We  start  with  the  equation  of 
motion  of  a  general  2nd  order  linear  system 


Kx  +  Cx  +  Me  =  / 

(5.4a) 

where 

x  =  Xe^^ 

(5.4b) 

II 

(5.4c) 

If  we  differentiate  Equation  (5.4b)  we  get  that 

X  =  JQXe^^  =  jQx 

(5.5) 

hence 

II 

(5^6) 

differentiating  Equation  (5.4b)  twice  we  get  that 

x=jQ(jQXe^)  =  JQx 

substituting  Equations  (5.6  )and  (5.7)  into  Equation  (5.4a)  we  get 

(5.7) 

+  Cx  +  MjQx  =  f 

(5.8) 

We  can  rewrite  equation  (5.8)  as 

( K  ^ 

—  +  C  +  MJQl  x  =  f 

(5.9) 

Following  reference  (5)  we  can  take  the  Laplace  transferm  of  equation  (5.9)  and 
thus  write  the  in^edance  of  the  linear  system  in  terms  of  the  complex  Laplace  parameter  s 
as 
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(5-10) 


Z{s)  =  Ms  +  C  +  — 
s 

In  reference  (5)  the  impedance  matrix  of  a  general  system  is  represented  as  an 
infinite  Laxirent  series  about  the  orgin  in  powers  of  the  con^lex  Laplace  parameter  s  as 

~  S  ~  '  '+A^s  +  ^  +  -+•  •  •+A_„  -^)  (5.11) 

n=-«  S  S 

After  defining  a  truncated  Laurent  series 

Z{s)  ^A,s  +  A,+A_,-  =  Ms  +  C+—  (5. 12) 

s  s 

which  approximates  the  impedance  matrix,  an  error  function  E(s),  and  a  cost  J  are  defined 
as 

E{s)  =  Z{s)-Z{s)  (5.13) 

and 

(5.14) 

where  W(s)  is  a  conplex  valued  weighting  fimction,  the  subscript  E  denotes  the  euclidian 
norm  of  a  matrix  and  the  integration  is  performed  over  a  prescribed  path  P  in  the  con^lex 
plane.  By  setting  the  partial  derivatives  of  J  with  respect  toM,Z,  andC  to  zero  the 
authors  of  reference  (5)  obtained  ejq)ressions  for  matricesji/,Z,  andC  which 
approximate  the  strfifiiess,  mass  and  damping  matrices  about  a  mode  of  the  system  The 
expression  for  M ,  Z ,  and  C  are  as  follows: 
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where 


K  = 


1 


ab~c^ 


a. 

c 

L 


V  ^1 

-  a  \—Zj{i^\W{i£t)\dO. 


(5-1 5c) 

a=  Jn2|jr(/Q)|c/Q 

(5-1 5d) 

(5-15e) 

c= 

(5-15f) 

2)  =  Z^(/Q)  +  /Z,(/Q) 

(5-15g) 

We  will  appfy  Equation  (5.15)  to  ^  as  defined  by  the  SST  and  obtain 

matrices  AAT  ,AM,  and  AC  which  will  serve  as  correction  matrices  of  the  FRF  of  the 
system  about  a  given  mode  of  the  system. 

As  an  example,  let  us  use  Equation  (5.15)  to  confute  the  (1,1)  element  of  the 
correction  matrix  AC  at  mode  1  (©1=25.21  Hz)  for  our  spatially  incomplete  system.  We 
take  the  weighting  fimction  to  be  W(in)=l/(in)  ,  the  path  P  is  taken  to  be  along  the 
imaginary  axis  fi-om  155.30  rads/sec  to  161.58  rads/sec. 


To  confute  this  line  integral  we  shall  use  the  trapezoid  rule  with  a  three  point 
fi’equency  sampling 

S  =  {155.30m<fc  /  sec,158.44racfe  /  sec,161.58rat&  /  sec}  (5.16) 

of  the  straight  line  pathP  (sampling  fi^quency  of  ti  rads/sec).  We  first  compute  the  value 
of  the  constant  a  of  Equation  (5.15d).  Using  vector  notation  we  have  that 

Q(H)  =  [I55.30rad:y  /  sec,158.44raffc  /  sec,161.58raffe  /  sec]  (5.17a) 

and 

W{a)  =  [0  -  0.0064y  sec/  rads,0  -  0.0063/ sec/  rads,0  -  0.0062 /sec/  rads]  (5.17b) 
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(5.18) 


hence  the  integrand  is  the  vector 

Q'(S)  .|)!r(S)|  =  [155.30,158.44,161.58] 

Usir^  the  trapezoid  rule  to  confute  the  integral,  we  have  that 

a=  j  =  +  ^  (519) 

155 30  2  V  2  2  y  /sec 

h=.  00000158  sec^/rads^  and  c—0. 0397  are  confuted  in  a  similar  manner. 

To  actually  compute  the  matrix  AC  ,we  use  the  SST  to  compute  the  matrices 

AZ(155.30  rads/sec),  AZ(155.44  rads/sec),  and  AZ(161.58  rads/sec)  for  the  points  of  the 

frequency  sarq)ling  of  the  straight  line  path  P.  Using  Equation  (2.26)  the  matrix 

AZ(1 55.30  rads/sec)  is  seen  to  be 


r(0.64E-9)  +  0.02y  (-0.04E  -  8)  -  0.04;  (0.03E-7)  +  0.01;  ' 
AZ(l55.44^^%^)=  (-0.4£-8)-0.04y  (-0.4£-9)-0.04y  (0.25£-7)-0.10y 

L(0.3E-8)  +  O.0iy  (025E-7)-0.10y  (-0.24E- 7)  +  0.15yJ 

(5.20) 

The  integral  in  Equation  (5.15b)  is  conqjuted  on  an  element  by  element  basis.  For 
(1>1)  dement  of  Equation  (5.15b),  if  we  collect  the  (1,1)  elements  of  the  matrices 
AZ(155.30  rads/sec),  AZ(1 55.44  rads/sec),  and  AZ(161.58  rads/sec)  and  use  the 
weighting  vector  given  in  Equation  (5.17b)  the  integrand  of  Equation  (5.15b)  is  the  vector 


AZ^‘’’^(H)|PF(S)|  =  [0.0064  0.040  -0.848][0.64  0.63  0.62](E-9)^^^se^ 


(5.21a) 


AZr(S)|r(S)|  =  [0.0041  0.2885  0.0036](E-9)^*-^se^^^  (5  2ib) 


Computing  the  integral  we  have 

AC(u)  =  -  J  Azr(S)|r(=«^=^f(^ 


AC(U)  =  =  oi??!  ^ 

(5.22a) 

AC(U)  =  0.0232£-6®-^  (5,22b) 

Performing  the  above  procedures  for  the  remaining  elements  of  Equation  (5.15b)  we  get 
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■  0.0232  -0.0443  0.0300 1 

AC  =  -0.0443  0.1278  -0.0914  (£-6)^^*^  (5.23a) 

0.0300  -0.0914  0.0645  J 

In  a  similar  manner  we  have 

■-0.0001  -0.0001  0.0002" 

AM  =  -0.0001  0.0015  -0.0020  /6/m  (5.23b) 

_  0.0002  -0.0020  0.0027  _ 

■-6.0476  4.5347  3.6757  1 

aZ=  4.5347  51.2196  -71.7061  (5.23c) 

_  3.6757  -71.7061  85.4417  J 

Figure  5-6  is  a  comparison  plot  of  experimental  and  uncorrected  analytic  FRFs 
versus  single  mode  integral  solutions  at  mode  1  using  a  frequency  sampling  of  3  points 
equally  spaced  over  bandwidths  of  25,  10,  and  1  Hz.  As  with  the  matrix  solutions,  the 
integral  solutions  are  sensitive  to  the  length  of  the  bandwidth  over  which  the  solutions  are 
computed  with  smaller  bandwidths  yielding  more  accurate  solutions.  Figure  5-7  is  a 
comparison  plot  of  ejqjerimental  and  imcorrected  analytic  FRFs  versus  single  mode 
integral  solutions  at  mode  1  over  a  1  Hz  bandwidth  using  frequency  samplings  of  3,10,50, 
and  200  points  from  the  bandwidth.  As  with  the  matrbc  solutions,  the  integral  solutions  are 
also  insensitive  to  the  sampling  frequency. 

Figures  5-8,  5-9,  and  5-10  are  comparison  plots  of  matrix  and  integral  solutions  at 
mode  1  computed  over  25,  10,  and  1  Hz  bandwidths  respectively,  using  three  point 
frequency  samplings.  These  three  plots  show  that  for  these  frequency  sampling,  the  matrix 
solutions  are  more  accurate  than  the  integral  solutions.  We  must  note  that  we  have  used 
the  trapezoid  rule  for  computation  of  the  integrals.  It  is  expected  that  better  accuracy 
would  be  achieved  from  the  integral  solutions  if  a  higher  order  integration  method  such  as 
Simpson's  rule  were  used. 
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Figure  5-  6  Spatially  incomplete  experimental  FRF  vs  single  mode  integral  solutions  at 
mode  1  using  3  point  frequency  samplings  of  25,  10,  and  1  Hz  bandwidths. 
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Figure  5-7  Spatially  incomplete  ejqjerimental  FRF  vs  single  mode  integral  solutions  at 
mode  1  using  3, 10,  50  ,  and  200  point  frequency  samplings  of  a  1  Hz  bandwidth. 


Omega  (Hz) 


Figure  5-  8  E5q)erimental  FRF  vs  single  mode  integral  and  matrix  solutions  at  mode  1 
using  a  25  Hz  bandwidth  with  a  3  point  frequency  sampling. 


Figure  5-  9  Ejqperimental  FRF  vs  single  mode  integral  and  matrix  solutions  at  mode  1 
using  a  10  Hz  bandwidth  with  a  3  point  frequency  sampling. 


60 


Figure  5-  10  Ejqjerimental  FRF  vs  single  mode  integral  and  matrix  solutions  at  mode  1 
using  a  1  Hz  bandwidth  with  a  3  point  frequency  sampling. 
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VI.  MULTIPLE  MODE  SOLUTIONS 

A.  MULTIPLE  MODE  MATRIX  SOLUTIONS 

We  wish  to  extend  the  results  of  Chapter  V  to  multiple  modes  of  the  system  imder 
consideration.  To  accomplish  this  we  shall  simply  extend  Equation  (5.1)  over  multiple 
frequency  samplings.  For  a  set  of  system  modes  {©i,  ()02,...,()i)„},  let  S((Bi)={Qii, 
Oi2,...,Oik}  be  a  frequency  sampling  of  a  bandwidth  [Qii,Qi„]  about  ©i  for  i=l,2,...,n.  We 
shall  apply  Equation  (2.30)  to  the  frequency  sampling 

“  =  (6-1) 

which  is  the  set  theoretic  union  of  the  sets  S(©i)  for  i=l,2,...,n.  If,  for  simplicy,  we  restrict 
ourselves  to  equally  sized  frequency  samplings  we  obtain  a  set  on  nk  equation  in  three 
unknowns 

■az,(q„)1  P  JQJ 

i  i  i  :  fAXf 

az,(q^)  =  I  -Qp  jap  am:  (6.2) 

:  :  i  i  [ACf_ 

[AKf 

One  can  solve  Equation  (6.2)  for  AA/f  as  we  have  done  in  Chapter  V  using  the 

.ACJ 

femihar  MATLAB  puedoinverse  function  but  better  results  are  obtained  if  we  weight  the 
equations  as  discussed  in  [Ref.  6].  A  good  weighting  is  to  assign  the  weight  ©i  to  those 
equations  associated  with  the  frequency  sampling  H(©i).  Alternately  we  could  assign  the 
weight  l/©i  to  those  equations  associated  with  the  frequency  sampling  S(©i).  The  ©i 
weighting  results  in  the  solution  being  more  accurate  at  the  lower  modes  while  the  l/©i 
weighting  gives  better  accuracy  at  the  higher  modes. 

We  note  that  the  solution  to  Equation  (6.2)  approximately  corrects  the  analytic 
FRF  in  the  sense  of  Equation  (5.1)  when  the  matrices  AZ(Qj)  as  given  by  Equation  (2.26), 
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are  computed  over  a  'good'  error  set,  Cerr.  In  general,  due  to  the  smearing  that  occurs 
during  reduction  of  the  full  analytic  system  to  the  reduced  analytic  system,  the  set  Cerr  is 
not  the  same  set  as  the  intersection  of  the  spatially  complete  error  set  and  the  A  set.  It  is  in 
general  a  larger  set  than  this  intersection  and  the  size  is  a  function  of  the  type  of  error,  i.e., 
mass,  damping,  or  stifi&iess  error  and  the  locations  of  the  errors.  Figure  6-1  is  a 
comparison  plot  of  experimental  and  uncorrected  analytic  FRFs  versus  the  multiple  mode 
matrix  solution  corrected  FRF  with  the  solution  having  been  computed  over  modes  1  and 
2  using  a  1  Hz  bandwidth  about  each  mode  with  a  3  point  frequency  sampling  over  the 
bandwidth  at  each  of  the  modes.  Figure  6-2  is  a  comparison  plot  using  a  matrix  solution 
computed  over  mode  1  through  4  again  using  a  1  Hz  bandwidth  about  each  mode  with  a  3 
point  frequency  sampling  over  the  bandwidth  at  each  of  the  modes. 
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H(9,9)  in/lbfaoglOof)  H(9,9)  in/lbf  GoglO  of) 


Figure  6- 1  Experimental  FRF  vs  multiple  mode  matrix  solutions  at  modes  1  and  2  using 
a  1  Hz  bandwidth  with  a  3  point  frequency  san^lings  at  each  mode. 
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Figure  6-  2  Experimental  FRF  vs  multiple  mode  matrix  solutions  at  modes  1  through  4 
using  a  1  Hz  bandwidth  with  a  3  point  frequency  samplings  at  each  mode. 
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B.  MULTIPLE  MODE  INTEGRAL  SOLUTIONS 


We  wish  to  extend  single  mode  integral  solutions  to  multiple  modes  of  the  system. 
As  with  the  multiple  mode  matrix  solutions  we  simply  take  as  the  path  of  integration  a 
path  which  is  the  set  theoretic  union  of  suitable  paths  at  each  of  the  modes  under 
consideration.  As  with  the  matrix  solutions  one  can  chose  the  weighting  function  W(s)  to 
be  1/(Q)  or  Q  to  achieve  lower  or  higher  mode  accuracy.  As  stated  in  the  previous  section 
it  is  required  that  a  'good'  error  set  is  known.  Figure  6-3  is  a  comparison  of  ejqserimental 
and  uncorrected  FRF  versus  multiple  mode  integral  solution  corrected  FRF  with  the 
solution  having  been  computed  over  modes  1  and  2  using  a  1  Hz  bandwidth  at  each  mode 
with  a  3  point  frequency  sampling  over  the  bandwidth  at  each  of  the  modes.  Figure  6-4  is 
a  comparison  using  a  solution  coir^uted  over  mode  1  through  4  again  using  a  1  Hz 
bandwidth  at  each  mode  with  a  3  point  frequency  sampling  over  the  bandwidth  at  each  of 
the  modes.  Figure  6-5  is  a  comparison  plot  of  the  multiple  mode  matrbc  and  integral 
solutions  over  modes  1  through  mode  4  using  a  1  Hz  bandwidth  at  each  mode  with  a  3 
point  frequency  sampling  over  the  bandwidth  at  each  of  the  modes. 
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H(9,9)  in/lbf(loglOof)  H(9,9)  in/lbf  OoglOof) 


Figure  6-  3  Experimental  FRF  vs  multiple  mode  integral  solutions  at  modes  1  and  2  using 
a  1  Hz  bandwidth  with  a  3  point  frequency  samplings  at  each  mode. 


67 


Figure  6-  4  Ejqjerimental  FRF  vs  multiple  mode  integral  solution  at  modes  1  through  4 


using  a  1  Hz  bandwidth  with  a  3  point  frequency  samplings  at  each  mode. 
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Figure  6-  5  E3q)erimeiital  FRF  vs  multiple  mode  matrix  and  integral  solutions  at  modes  1 
and  2  using  a  1  Hz  bandwidth  with  a  3  point  frequency  samplings  at  each  mode. 
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C.  SINGLE  POINT  MULTIPLE  MODE  SOLUTIONS 


The  results  of  sections  A  and  B  can  be  performed  using  partitions  consisting  of  a 
single  point.  Figure  6-6  shows  the  results  of  performing  a  multiple  mode  matrix  solution 
over  the  first  three  modes  of  our  spatially  incomplete  beam  using  a  single  point  partition  at 
each  of  the  modes  under  consideration.  Figure  6-7  shows  the  results  of  performing  a 
multiple  point  integral  solution  over  the  first  three  modes  of  our  spatially  incomplete  beam 
using  a  single  point  partition  at  each  of  the  modes  under  consideration.  Figure  6-8 
compares  the  multiple  mode  matrix  and  integral  solution  computed  over  the  first  three 
modes  of  our  spatially  incomplete  beam  with  the  analytic  FRF.  In  figure  6-9  we  compare  a 
multiple  mode  matrix  solution  computed  over  the  first  three  modes  of  a  spatially  coinplete 
beam  using  a  single  point  partition  at  each  of  the  modes  with  the  beam's  spatially  complete 
analytic  FRF.  In  figure  6-10  we  compare  a  multiple  mode  integral  solution  computed  over 
the  first  three  modes  of  a  spatially  complete  beam  using  a  single  point  partition  at  each  of 
the  modes  with  the  beam's  spatially  complete  analytic  FRF.  Figure  6-11  compares  multiple 
mode  matrix  and  intgral  solution  computed  over  the  first  three  modes  of  a  spatially 
complete  beam  with  the  beam's  spatially  complete  analytic  FRF.  Figures  6-9  and  6-10 
show  that  the  single  point  multiple  mode  solutions  are  exact. 
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H(9,9)  in/IbfaoglOof)  H(9,9)  in/lbf  aoglOof) 


Figure  6-  6  E3q)erimental  FRF  vs  multiple  mode  matrix  solution  computed  at  modes  1 
through  3  using  a  1  Hz  bandwidth  with  a  single  point  frequency  samplings  at  each  mode 
for  a  spatially  incomplete  beam. 


i 


I 
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H(9,9)  in/lbf(loglOof)  H(9,9)  in/lbf  OoglOof) 


Figure  6-  7  Experimental  FRF  vs  multiple  mode  integral  solution  computed  at  modes  1 
through  3  using  a  1  Hz  bandwidth  with  a  single  point  frequency  samplings  at  each  mode 
for  a  spatially  incomplete  beam. 
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Figure  6-  8  Experimental  FRF  vs  multiple  mode  matrix  and  integral  solutions  computed  at 
modes  1  through  3  using  a  1  Hz  bandwidth  with  single  point  frequency  samplings  at  each 
mode  for  a  spatially  incomplete  beam. 
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H(6,6)  in/lbf  OoglO  of)  H(6,6)  in/Ibf  aoglO  of) 


Figure  6-  9  Experimental  FRF  vs  multiple  mode  matrix  solution  computed  at  modes  1 
through  3  using  a  1  Hz  bandwidth  with  a  single  point  frequency  samplings  at  each  mode 
for  a  spatially  complete  beam. 
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H(6,6)  in/lbf  (loglO  of)  H(6,6)  in/lbf  (loglO  of) 


Figure  6-  10  Experimental  FRF  vs  multiple  mode  integral  solution  computed  at  modes  1 
through  3  using  a  1  Hz  bandwidth  with  a  single  point  frequency  samplings  at  each  mode 
for  a  spatially  complete  beam. 
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H(6,6)  inAbfaoglOof) 


Figure  6- 11  Experimental  FRF  vs  multiple  mode  matrix  and  integral  solution  computed 
at  modes  1  through  3  using  a  1  Hz  bandwidth  with  a  single  point  frequency  samplings  at 
each  mode  for  a  spatially  complete  beam. 
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H(9,9)  in/lbf  (logic  of) 


Figure  6- 12  Experimental  FRF  vs  multiple  mode  matrix  solutions  computed  over  1,  2,  3, 
and  4  modes  using  a  1  Hz  bandwidth  with  single  point  frequency  samplings  at  each  mode 
for  a  spatially  complete  beam. 
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vn.  CONCLUSIONS  /  RECOMMENDATIONS 


A.  SUMMARY 

Frequency  Doinaiii  Structural  Identification  using  the  structural  s5mthesis 
transformation  was  performed  on  a  simulated  e3q)erimental  free-fi-ee  beam  model  The 

identification  was  performed  for  both  the  spatially  conqilete  and  the  spatially  inconq)lete 
case. 

•  SST  based  structural  identification  provides  an  exact  solution  for  the 
identification  of  FE  modeling  errors  given  spatially  complete  data. 

•  For  spatially  inconqilete  systems  SST  based  structural  identification  provides 
a  firequency  dependent  solution  which  is  not  suitable  for  finite  element  modeling  error 
correcting. 

•  Using  both  matrix  and  integral  based  techniques  the  SST  can  provide  singift 
mode  solutions  which  are  fi-equency  independent  correction  matrices  AK,  AM,  AC 
which  approximate^  corrects  the  reduced  FE  model  within  a  frequency  bandwidth  of 
the  ejqjerimental  system  mode  under  consideration. 

B.  CONCLUSIONS 

This  thesis  has  clearly  demonstrated  that  for  spatially  incon^lete  systems,  single 
mode  frequency  independent  solutions  can  be  found  which  correct  the  reduced  finite 
element  model  in  a  neighborhood  of  a  given  mode  of  the  experimental  system.  It  has  also 
been  shown  that  the  concept  of  a  single  mode  solution  can  be  expanded  to  that  of  multiple 
mode  solutions  which  are  frequency  independent  correction  matrices  which  approximately 
corrects  the  reduced  FE  model  throughout  a  frequency  bandwidth  which  includes  more 
than  1  mode  of  the  e3q)erimental  system. 

C.  RECOMMENDATION 

The  purpose  of  this  thesis  was  to  find  frequency  independent  multiple  mode 
solutions  which  coxild  be  used  to  approximately  correct  reduced  finite  element  models. 
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Although  satisfectory  results  were  obtained,  investigation  is  still  required  in  the  following 
areas: 


•Determine  the  relationship  between  errors  in  a  full  FE  model  and  those  of  an 
associated  reduced  FE  model. 

•Determine  a  method  to  'pullback*  reduced  FE  model  correction  to  full  order 
FE  model  corrections. 

•Investigate  use  of  the  integral  formulation  of  reference  (5)  in  analysis  of  the  O- 
set  system. 
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APPENDIX 


The  following  is  a  brief  description  of  MATLAB  routines  employed  in  this  thesis: 

•SST.M  -  Generates  MAT  files  containing  spatially  complete  or  spatially  incomplete 
experimental  FRF,  analytic  FRF,  localization  matrix,  and  SST  solutions. 

•PLOTSST.M  -  Uses  files  produced  by  SST.M  to  generate  plots  used  in  chapters  3  and  4. 

•CHAP3.M  -  Sets  parameters  and  calls  SST.M  and  PLOTSST.M  to  generate  figures  for 
chapter  3. 

•CHAP4.M  -  Sets  parameters  and  calls  SST  M  and  PLOTSST.M  to  generate  figures  for 
chapter  4. 

•CHAP5.M  -  Generate  figures  for  chapter  5. 

•CHAP6.M  -  Generate  figures  for  chapter  6. 

•SETUP.M  -  Generates  MAT  files  for  FE  models. 

•BEAMMDL.M  -  Setup  FE  models. 

•FSTATIC.M  -  Perform  static  reduction  of  mass  and  stiffness  matrices. 

•FIRS_TAM.M  -  Perform  IRS  reduction  of  mass  and  stiffhess  matrices. 

•FREQMODE.M  -  Returns  FE  model  fi-equencies. 

•NDX3D.M  -  Indexes  a  series  of  2-D  matrices  into  a  single  3-D  matrix. 

•INTSUB.M  -  Decomposes  impedance  matrix  into  mass,  stiffhess,  and  damping  matrices  using 
integral  techniques. 

•MYTRAP.M  -  Computes  integrals  using  trapezoidal  method 
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SST.M 


%%%%%General  setup%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%%%%%%%0/o0/o0/o%%0/o0/o0/o0^0^0/^0/^0/^0/^0/^0^^^^^ 

clear,clc  ;%%clear  workspace. 

5  closeall  ;%close  aiiy  open  figure  windows 

j=sqrt(-l)  ; 

load  sstconf  ;%load  A-set  and  O-set 
aset(oset)=[]; 

1 0  aset_size=iengtli(aset); 

oset_size=lengtli(oset); 
load  beamdata 
if  lengtli(oset)==  0 
complete=l; 

15  else 

coinplete=0; 

end 


%%setup  plot  Iabels%%%%%%%%%%%%o/o%o/o%o/o0/o0/o0^o/^o/^o/^o^o^o/^o/o0/^^^^^^ 
20  if  (oset=[])&(posofMasserr(  1  >^0)&(posofStifferr(  1  >-=0) 

casenanie=(’Spatially  Complete  Mass  &  Stiffiiess  error’); 
elseif  (oset==[])&(posofMasseiT(l  >^0) 
casenanie=(’SpatialIy  Complete  Mass  error'); 
elseif  (oset==[])&(posofStifrerT(l  >-=0) 

25  casename=(’Spatially  Complete  Stifliiess  error*); 
elseif  (posofMasserrf  1  )-=0)&(posofStitTerr(  1  >-^0) 
casename=('Spatially  Incomplete  Mass  &  Stiffness  error'); 
elseif  (posofMasserrf  1  )-=0) 
casename=('Spatially  Incomplete  Mass  erroi^); 

30  elseif  (posofStifferr(  1  >-=0) 

casename=('Spatially  LiComplete  Stiffness  error'); 
end 

%%get  reduced  A  set  &  O  set  frequencies%%%%%%%%%%%%% 

35 

kexto=k_anal(oset,oset);  %%  stiffiiess 
mexto=m__anal(oset,oset);  %%  and  mass  matrices 
if  complete 
kstat=k_anal; 

40  mstat=m_anal; 
else 

if  static  “  0 

[kstat,mstat]=fstatic(k_anal,m__anal,oset,aset);%%  get  reduced  K  &  M 
elseif  static  ==  1 

45  [kstat,mstat]=firs_tam(k_anal,m_anal,oset,aset);%%  get  reduced  K  &  M 

else 

kslat=k_anal(aset,aset); 

mstat=m_anal(aset,aset); 

end 

50  end 

cstat=sqrt(- 1  )*struc_damping.  *kstat; 

[u,lamb^a,c]=freqmode(k_anal,m_anal)  ;%get  freqs 
[u,lambdared,c]=freqmode(kstaMnstat) ;% 

[u,lambdax,c]=freqmode(k_exp,m_exp)  ;% 

55  [u,lambdaexto,c]=ft^eqmode(kexto,mexto);% 

omegaa=sqrt(lambdaa); 
omegax=sqrt(lambdax); 
omegared=sqrt(lambdared); 
omegaex1o=sqrt(Iambdaexto); 

60  for  i=l  :4  %  tliis  to  delete  fixed  body  modes 
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SST.M 


65 


70 


75 


ifreal(omegax(l))<  10^(-3) 

omegax=oinegax(2:length(omegax)); 

end 

ifreaI(omegaa(l))<  10^(-3) 
omegaa=omegaa(2:lengai(omegaa)); 
end 

ifcomplete==0 

if  real(omegared(  1 ))  <  1 0^(-3) 
omegared=omegared(2:Iengtli(omegared)); 

end 

if  reaI(omegaexto(  1 ))  <  10^(-3) 
omegaexlo=omegaexlo(2:lengai(omegaexlo)); 

end 

end 

end 
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00 


15 
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%%%%%fmd  true  mass  and  stiffiiess 
true_stiflhess=k_exp-k_anal; 
stifBiess_cset=find(diag(tnie_stifrness)); 
true_niass=m_exp-in_anal; 
mass_cset=find(diag(true_mass)); 
true_damping=c_exp-c__anal; 
dainping_cset=find(diag(true_dainping)); 
if  lengtli(mass_cset)  >  0 
for  i=l:lengtli(mass_cset) 

x=fmd(stiflhess_cset==mass_cset(i)); 

if'-<lengtli(x)>0) 

stifBiess_cset=[sdffness_cset;  mass_cset(i)]; 
end 
end 
end 

if  lengtli(damping_cset)  >  0 
for  i=l:lengtli(damping_cset) 

x=find(stijf]&iess_cset==^daniping_cset(i)); 

if-<lengtli(x)>0) 

stifliiess__cset=[stifliiess_cset;danipmg_csel(i)]; 

end 

end 

end 

true_cset=sort(stiffiiess_cset) 

save  true_errs  true_cset  tnie_mass  true_stifihess  true_daiiiping 
clear  true  jxiass  true__stiffiiess  tme_daniping 

clear  area  eeii  force  g  gc  gcx  gk  gkx  gm  griLX  goblc  goblcx  goblk  goblkx 
clear  goblm  goblnix  ke  lambdaa  lainbdaexto  lambdared  lambdax  lumpdamp 
clear  limipinass  lunipspring  me  pho  posofDamperr  posofMasserr  posolStifferr 
clear  stiffiiess  cset 


%%%setup  frequency  rmige  for  aiialysis%%«/oO/.«/„o/oO/oO/„o/,o/„o/„.^^^^^^^ 

freqtop=omegax(luglmiode)+.  1  *omegax(lughmode);  o/o/o/o/o/o/o/o 

freqbottom=omegax(lowmode)-.  1  *omegax(lowmode); 

freqbottom=10*2*pi;  %startat  lOliz 

freqtop=1200*2*pi;  %end  at  1200  liz 

%%nmnber  of  points  to  plot 

nmnpoints=400;%%o/o%o/oSET  NUMBER  OF  POINTS  TO  USE  FOR  THE  SIMULATION 
imeness=nuinpoints; 


w=freqbottom:(freqtop-freqbottom)/(numpoints-l);freqtop; 

%%l^CateCset%%%%%o/oO/„o/„o/„o/„o/„o/„o/,o/,o/„o/„o/„o/„o/„o/„«/„o/„o/,«/^^^^^ 
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SST.M 


keepgomg=’n’ 
wlreq=omegax(lo\viiiode) ; 
iiextinode=lowiiiode; 
lol=l  ; 

1 25  defaultstepsize=.  1 ; 

while  (keepgoing~'n’) 

z_aiial_red=kstat+j  ♦  wfreq*cstat-wfreq^2  *mstat ; 

h_aiial_red=inv(z_anal_red) 

z_exp=k_exp+j  *wfreq*c_exp-wlreq'^2*in__exp  ; 

130  h_exp=mv(z__exp)  ; 

h__exp=h_exp(aset,aset)  ; 

L=z__aiial__red*(h__anal_red-h_exp)*z__anal_red  ; 

ci=fhid(abs(diag(L))>tol) 

lemp_cset_size=lengtli(ci); 

1 35  if  leiigtli(ci)  >=  1 

cset_size=lengtli(ci) 

figure(l) 

plot(aset,  abs(diag(L)))  ; 
hold  on 

1 40  plot(aset(ci),zeros(  1  ,lengtli(ci))/x’) 

tallest=max(abs(diag(L)))  ; 

title([TOL  =  ',nuni2str(tol),'  Cset  size  =  ’,int2str(cset_size),'  Oinega  =  ',iuun2stri(wtreq/(2*pi)),'  Hz’]) 

ylabel(’ibf7in') 

if  complete 

145  xlabel(’DOF) 

else 

xlabel(’ASETDOF) 

end 

hold  off 

150  k=meuu(  '  Choose  an  action 

Increase  Tolerance 

Decrease  Tolerance 
Increase  frequency 
Decrease  frequency 

155  ’  Previous  mode 

’  Next  Mode 
'Set  Print  switch  on’,... 

’  hicrease  stepsize 
’  Decrease  stepsize 

160  ’Go  ’); 

ifk==l 

stepsize=defaultstepsize; 

ctr=0; 

wliile  lengtli(ci)  >=  temp_cset__size 
165  ctr=ctrM; 

ifctr==10 

stepsize=stepsize*  1 0; 
ctr=0 ; 
end 

170  toI=min(abs(tol+stepsize*(tol)), tallest) ; 

ci=fuid(abs(diag(L))>tol)  ; 

iflengtli(ci)>=  1 
cset_size=length(ci)  ; 

figure(l) 

175  plot(aset,  abs(diag(L)))  ; 

hold  on 

plot(aset(ci),zeros(  1  ,lengtli(ci)),'x’) 
taIlest=max(abs(diag(L)))  ; 

tille([TOL  =  ',nuni2sfr(tol),’  Cset  size  =  ',int2str(cset_size),'  Omega  =  ',nuni2str(\\freq/(2*pi)),’  Hz’]) 
180  ylabelCIbf/in') 
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185 


190 


195 


200 


205 

Hz’]) 

210 


1215 


(20 


(25 


!30 


35 


■40 


if  complete 

xlabel(T)OF') 

else 

xlabelCASETDOF) 

end 

hold  off 
end 
end 

elseifk==2 

stepsize=defaultstepsize; 

ctr=0; 

while  length(ci)  <=  lemp_cset_size 
ctr=clrf-l; 
ifctr==10; 
stepsize=stepsize*2; 
cti^O; 
end 

tol=niax(abs(tol-stepsize*(toI)),.0000 1 ) ; 
ci=fmd(abs(diag(L))>tol)  ; 
iflengtli(ci)>=  1 

cset_size=lengtli(ci)  ; 

figure(l) 

plot(aset,  abs(diag(L)))  ; 
hold  on 

plot(aset(ci),zeros(  1  ,lengtli(ci)),’x’) 
lallest=max(abs(diag(L)))  ; 

title([TOL  =  ',num2str(tol);  Cset  size  =  ’,int2slr(cset_size);  Omega  =  ’,num2str(wfreq/(2*pi)); 

ylabel(’lbf/in’) 
if  complete 
xlabelCDOF') 
else 

xlabelCASET  DOF’) 
end 

hold  off 
end 
end 

elseifk==3 

while  lengtli(ci)  ==  temp_cset_size 

wfreq=mm(wfreq+.0000 1  *(lreqtop-wfreq),freqtop); 
z_anal_red=kstat+j*w(comit)*cstat-wlieq^2*mstat ; 
h_anal_red=inv(z_anal_red) 
z_exp=k_exp+j*wfreq*c_exp-vvfreq^2*m_exp  ; 
h_exp=inv(z_exp) 
h_exp=h_exp(aset,aset)  ; 

L=z_anal_red*(h__anal_red-h_exp)*z_anal_red  ; 

tallest=max(abs(diag(L)))  ; 
tol=min(tol+.001*tallest,tallest) ; 
ci=fmd(abs(diag(L))>tol)  ; 
cset_size=lengtli(ci)  ; 

figiire(l) 

plot(aset,  abs(diag(L)))  ; 
hold  on 

plot(aset(ci)^eros(  1  ,lengtli(ci)),’x') 

title([TOL  =  ’,num2str(tol);  Cset  size  =  ’,int2str(cset_size);  Omega  =  ',num2str(wfreq/(2*pi));  Hz’]) 

ylabel(’lbf7m') 

if  complete 

xlabel('IX)F') 

else 

XlabelCASET  DOF') 
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'  end 
hold  off 
end 

elseifk==4 

245  wfreq=iiiax(wfreq-,05  *(  wfreq-freqbottom),freqbottoui); 

elseifk==5 

nextmode=niax(lowniode,nexUuode-l ); 
wfreq=omegax(nextinode); 
elseifk==6 

250  nexUnode=miii(liiglmiode,nextniode+ 1 ); 

wfreq=omegax(nextinode); 

elseifk==7 

pswitcli=y  ; 

elseifk==8 

255  defaultstepsize=defaultstepsize/10; 

elseifk==9  - 

defaultstepsize=defaultstepsi2e*  1 0; 


else 

keepgoing=y 

260 

cset=aset(ci); 

cset_rel=ci 

end 

else 

toI=tol*.9; 

265 

end 

end 

close(l) 

save  extrtset  cset  cset_rei 
%%% 

%%%%%%%%%%%%%%%%%%FORCE  THE  CSET%%%%%%%%%%%%%%%%%%%%%%o/oO/oO/o%o/qO^%o/^,^^ 
%cset=[7  9  11  13]  %%%%%%%%THE  INCOMPLETE  CASE  CSET%%%%%%%%%%%% 

%cset_rel=[4  5  6  7]%%%%%%%%%%%%%%%%%%%%%%%o/o%%o^%%o^%%o^o^y^,o^y^JO^%y^,o/^o^o^^ 

%%% 

%cset=[4  5  6  7  8  9  10  11]%%%%%%%%%%%%%THE  COMPLETE  CASE  CSET%%% 

%cset__rel=[4  5 6  7 8 9  10  ll]%%%%%%%%%%%%%%%%%%%%%%%%%%%%%o/oO/oO/oO/^o/^ 
%%%%%%%%%%%%%%%%%%%%%%yo%%%%%%%%%%%%%%%%%%%%o/o%o/o%%%%o/o%o^o^o/o%%%%% 
280  %%% 

cset_size=length(cset); 
save  L  L 
clear  L 

clear  kexta  kexto  niexta  mexto 
285  pack 

skyfull=zeros(nmndof,nuindof); 
count=l  ; 

for  indexl=l:muiidof 
for  index2=indexl:nunidof 
290  skyfull(mdexl,index2)=count; 

skyfull(index2,indexl  )=count; 
count=count+l; 
end 
end 

295  skyred=zeros(aset_si2e,aset_size); 
count=l  ; 
for  indexl=l:aset_size 
for  index2=mdexl:aset_si2e 

skyred(indexl  ,index2)=count; 

300  skyred(mdex2,indexl)=count; 
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coimt=couiit+l; 

end 

end 

skycset=zeros(cset_size,cset_size); 

305  count=l  ; 

for  index  1=1  :cset_size 
for  iiidex2=indexl:cset_size 

skycset(indexl  ,iudex2)=couiit; 
skycsel(index2,mdexl  )=coimt; 

310  count=count+l; 

end 
end 

fiiiniolder=zeros(  1  ,(nuindojf^(nunidof+ 1  ))/2); 

315  redJioldei=zeros(  1  ,(aset_size*(aset_size+l  ))/2); 

cset Jiolder=zeros(  1  ,(cset_size*(cset_size+ 1  ))/2); 

if  meters 

\vaitbarJiandle=waitbar(0, 'Computing  Experimental  FRF'); 

320  else 

dispCGetting  experimental  FRF) 
end 

H_EXP=Q; 

for  count=  1  :numpoints 
325  if  meters 

waitbar(count/numpoints); 

end 

z_exp=k_exp+j*w(count)*c_exp-w(count)^2*m_exp; 

h_exp=mv(z_exp); 

330  h_exp=h_exp(aset,aset); 

skyindex=l; 
for  index  1 = 1  laset^size 

for  mdex2=indexl:aset_size 
red_holder(skyindex)=li_exp(indexUindex2); 

335  skyindex=skyindex+l ; 

end 
end 

H_EXP=[H^EXP  redjiolder^]; 
end 

340  if  meters 

close(waitbar_handle) 

end 

saveH^EXPH_EXP 
clear  H__EXP 
J45  pack 

if  meters 

\vaitbarJiandle=waitbar(0, 'Computing  Experimental  Impedence'); 
else 

50  dispCGetting  experimental  Impedence') 

end 

Z^EXP=[]; 

for  count=l:numpomts 
if  meters 

55  waitbar(couiit/numpoints); 

end 

z_exp=k_exp+j*w(count)*c_exp“Vv(count)^2*ni_exp; 

h_exp=mv(z_exp); 

h_exp=h_exp(aset,aset); 

30  z_exp=inv(h_exp); 
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365 


370 


375 


skyindex=l ; 

for  index  1 = 1  :aset_size 

for  iiidex2=indexl:aset_size 
red JioIder(skyindex)=z_exp(index  1  ,index2); 
skyindex=skymdex+ 1 ; 
end 
end 

Z_EXP=[Z^EXP  redjiolder']; 
end 

if  meters 

close(\vaitbar_handIe) 

end 

saveZ_EXPZ_^EXP 
clear  Z_EXP 
pack 


if  meters 

waitbar  JiandIe=waitbar(0;Computing  Reduced  Analytic  FRF'); 
else 

380  dispCGetting  Reduced  FRF’) 
end 

H_ANAL_RED=[]; 
for  count=l:numpoints 
if  meters 

waitbar(count/numpoints); 

end 

z_anal_red=kstat+j*w(count)*cstat-w(coimt)^2*mstat; 

li_anal_red=inv(z_anal_red); 

skyindex=l  ; 

390  for  index  1=1  :aset_size 

for  iudex2=indexl:aset__size 
redJiolder(skyindex)=li_anal_red(indexl,index2); 
skyindex=skyindex+ 1 ; 
end 

395  end 

H_ANAL^REI>[H^ANAL  RED  red  holder^]- 
end 

if  meters 

close(\vaitbar_handle) 

400  end 

save  H_ANAL_R  H_ANAL_RED 

clear  Z_ANAL_RED  H_ANAL_RED  kexta  kexto  mexla  mexto 
pack 


405  %%  compute  and  save  L  matrix  diagonal  as  a  function  of  frequency%%%%%%%% 

if  meters 

waitbarJiandle=waitbar(0;Coniputing  L  Diagonals'); 
else 

dispCGetting  L  Matrix') 

410  end 

L^DIAGS=[]; 
for  count=l:numpomts 
if  meters 

vvaitbar(count/numpoints); 

415  end 

z_anal=k_anal+j*w(count)*c_anal-w(count)^2*m_anal; 
z_ana!_red=kstat+j  *  w(count)*cstat-’W(count)'^2  *mstat; 
h_aual_red=inv(z_anal_red); 
z_exp=k__exp+j*c_exp-\v(count)"'2*m_exp; 

420  h_exp=inv(z_exp); 
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h_exp=h_exp(aset,aset); 

L=z__anal__red*(h_janal_red~li_exp)*z_aiial_red; 

tempJ_diags=diag(L); 

L_DIAGS=[L_DIAGS  temp__l_diags(:)]; 

425  end 

if  meters 

close(waitbar_handle) 

end 

save  L^DIAGS  L^DIAGS 
430  clear  L_DIAGS  z_anal  temp_I_diags  L 

clear  z__exp  z_anal_red  u  v  h_exp  h_anal_red 
pack 


435 


440 


445 


450 


455 


1460 


P65 


\70 


75 


!0 


if  meters 

\vaitbarJiandle=waitbai(0;Coinputing  Z  Analytic'); 
else 

disp('Getting  Z  anal') 
end 

Z_ANAL=[]; 
for  count=  1  rnumpoints 
if  meters 

waitbar(count/numpoints); 

end 

z_aiid=k_aiialHj*\v(count)*c_mial-vi<couiit)^2*m_aiial; 

skyindex=l  ; 

for  indexl=l:numdof 

for  mdex2=indexl:numdof 
fullJiolder(skyindex)=z_aiial(indexl,iiidex2); 
skyindex=skyindex+l ; 
end 
end 

Z^A]SfAL=[Z^ANAL  fulljioldef]; 
end 

if  meters 

close(waitbarJiandle) 

end 

save  Z_ANAL  Z_ANAL 
clear  Z__ANAL  z_anal 
pack 
if  meters 

waitbar__handle=vvaitbar(0;Computing  Reduced  Z  Anal>1ic'); 
else 

disp('Getting  Reduced  Z  anal') 
end 

Z^ANAL_RED=[]; 
for  count=l  :mmipoints 
if  meters 

waitbar(count/numpoints); 

end 

z_anai_red=kstat+j*w(count)*cstat-\v(count)^2*mstat; 

skyindex=l  ; 

for  indexl  ==  1  :aset_size 

for  index2=indexl:aset_size 
redJiolder(skyiiidex)=z_anal_red(indexl  ,index2); 
skyindex=skyiudex+l ; 
end 
end 

Z_ANAL_^RED=[Z^ANAL_RED  redjiolder']; 
end 

if  meters 
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close(waitbar_haiidle) 

end 

save  Z^ANAL^R  Z_ANAL_^RED 
clear  Z_ANAL__RED  z_anal_red 
485  pack 

%%%COmputeDZ%%%%%%%%%%%%%%%%%%%%%%o/oO/o%o/oO/oO/,o/„o/,o/„o/„o/„o/,o/„o/,o/oO/„o/,o/„o/,o/„o/„o/,o/„o/,^^^^^ 
if  meters 

waitbar Jiandle=waitbar(0,'Computmg  DZ. . 

490  else 

disp(’Getting  DZ') 
end 
DZ=[]; 

for  count=l  :numpoints 
495  if  meters 

waitbar(count/nmnpoints); 

end 

z_anal_red=kstat+j  *  w(count)*cstat-\v(count)''2  *mstat; 
h_anal_red=inv(z_anal_red); 

500  z__exp=4:_expf j  *\v(count)*c_exp-\v(count)^2  ♦m__exp; 

h_exp=mv(z_exp); 
h__exp=h_exp(aset,aset); 
hacc=h_anal_red(cset_rel,cset_rel); 
lixcc=h_exp(cset_rel,cset_rel); 

505  dz=inv(mv(mv(hacc)*(liacc-hxcc)*inv(hacc)>-hacc); 

skyiudex=l  ; 
for  index  1 = 1  :cset_size 

for  index2=indexl:cset_size 
cset Jioider(skymdex)=dz(index  I  ,index2 ); 

510  skyindex=sl^index+l ; 

end 
end 

DZ=[DZ  csetjiolder*]; 
end 

515  if  meters 

close(waitbarjiandle) 

end 

save  DZ  DZ 

clear  hacc  hxcc  dz  h_exp  z_anal_red  h_anal_red  z_exp 
520  pack 

%%%decompose  DZ  into  DM,  DK,  &  DC%%%%%%%%%%®^%%%«^%o^%o/oO/oO/o0^o/„o^o^o/^,o/„o^o/,o^o^o/„o/,%%«^o^o/„ 
if  meters 

waitbar Jiandle=\vaitbar(0,T)econiposing  DZ. . 

525  else 

dispCDecomposing  DZ') 
end 
DK=[]; 

DM=D; 

530  DC=n; 

dz  1  =zeros(cset_size,cset_size); 
dz2=zeros(cset_size,cset_size); 
dz3=zeros(cset_size,cset__size); 
for  i=l  :numpoints-2 
535  if  meters 

waitbar(i/(nmiipoints-2)) 

end 

dztemp=DZ(ndx3d([l  cset_size*(csel_sizei- 1  )/2  numpoints],  1 , 1  :cset_size*(cset_si2ef  1  )/2,i)); 
sk7index=i; 

540  for  indexl=l:cset  size 
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for  index2=mdexl:cset_size 
dzl  (index  1  ,index2  )=dztemp(skyindex); 
dzl  (index2,indexl  )=dztenip(skyindex); 
skyindex=skyindex+l ; 
end 
end 

d2temp=DZ(ndx3d([l  cset_size*(cset_sizef  1  )/2  niunpoinls],  1 , 1  :cset_size*(cset_sizef  I )/2,i+l )); 

skyindex=l; 

for  index  1 = 1  :cset_size 

for  index2=indexl:cset_size 
dz2(index  1  ,index2)=dztemp(skyindex); 
dz2(index2 , index  1  )=dztemp(skyindex); 
sky  index=skyindex+ 1 ; 
end 
end 

dztemp=DZ(ndx3d([l  cset_size*(cset_size+ 1  )/2  numpoints],  1 , 1  :cset_size*(cset_sizefl  )/2,i+2)); 

skyindex=l; 

for  indexl=l:cset_size 

for  index2=indexl:cset_si2e 
d23(indexl  ,index2)=dztemp(skyindex); 
dz3(index2,indexl  )=dztenip(skyindex); 
skyindex=skyindex+ 1 ; 
end 
end 

temp=[eye(cset_size)  -w(i)'^2*eye(cset_size)  +j*w(i)*eye(cset>ze); 
eye(cset_size)  ’•w(i+l)‘^2*eye(cset_,size)  +j*w(i+l)*eye(cset_size); 

eye(cset__size)-w(i+2)^2*eye(cset_size)+j^\v(i+2)*eye(cset_^size)]\[dzl;dz2;dz3]; 

ktemp=temp(  1  :cset_size,:); 
mtemp=temp(cset__size+l;2*cset__3ize,:); 
ctemp=temp(2*cset_si2e^’l  :3*cset_size,:); 

skymdex=l  ; 
for  index  1 = 1  :cset_size 

for  index2=indexl  :cset_size 
cset_holder(skyindex)=ktenip(indexl,index2); 
skyindex=skyindex+l ; 
end 
end 

DK=pK  cset__holder(:)]; 

skymdex=l  ; 
for  indexl = 1  :cset_size 

for  index2= indexl  :cset_size 
cset Jiolder(skyiiidex)=mtemp(index  1  ,index2 ); 
skyindex=skyindex+l ; 
end 
end 

DM=pM  csetJiolder(:)]; 

skyindex=l  ; 
for  indexl =1  :cset_size 

for  index2=index  1  :cset_size 
csetJiolder(skyindex)=ctemp(indexUiiidex2); 
sky  index=sky  index+ 1 ; 
end 
end 

DC=[DC  cset_holder(:)]; 
end 

if  meters 

close(vvaitbar_handle) 
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end  ‘ 
save  DM  DM 
save  DK  DK 
save  DC  DC 

605  clear  DC  DK  DM 

save  exp  c_exp  k_exp  m_exp 
clear  c_exp  k_exp  m_exp 
save  stat  kstat  mstat  cstat 
clear  kstat  mstat  cstat 

610  save  anal  k_anal  m_anal  c__anal 

clear  k_anal  m__anal  c__anal 
pack 


615 


620 


625 


630 


635 


640 


645 


650 


655 


660 


if  complete 
load  Z_ANAL 
Z=Z^ANAL; 
clear  Z^ANAL 
else 

loadZ_ANAL__R 
Z=Z_ANAL_RED; 
clear  Z__ANAL_RED 
end 

CORRHA=G; 

tempza=zeros(aset_size,aset__size); 
tempz=zeros(cset_size,cset_size); 
if  meters 

waitbar jiandle=waitbar(0,'lnstalling  DZ. . 
else 

disp(’Installing  DZ’) 
end 

for  i-l:numpoints 
if  meters 

waitbar(i/nmupoints) 

end 

ztemp=Z(ndx3d([l  aset_size*(aset_size+l)/2  (niuiipomts-(i-l))],l,l:aset  size*(aset  size+lV2  D)- 
skyindex=l;  ”  wy, 

for  indexl=l  :aset_size 

for  index2=indexl  :aset_size 
tempza(indexl  ,index2)=ztemp(skymdex); 
tempza(iiidex2,indexl  )=ztemp(skyindex); 
skyindex=skyindex+l ; 
end 
end 

ztemp-DZ(ndx3d([l  cset_size*(cset_si2e+iy2  (numpomts-<i-l ))],!, Ixset  size*(cset  size+lV2  nv 
skyindex=l;  ""  “  ^ 

for  index  1 = 1  :cset_size 

for  index2=indexl:cset_size 
tempz(index  1  ,index2  )=ztemp(skymdex); 
tempz(index2,indexl  )=ztemp(skyindex); 
skyindex=skyindex+l ; 
end 
end 

tenipza(cset__rel,cset_rel)=tempza(cset  rel,cset  relVf-tempz- 
DZ(:,1)=Q;  ~ 

tempha=mv(tempza); 

CORRHA=[CORRHA  tempha(:)]; 
end 

if  meters 

close(waitbar__handle) 
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end 

clear  DZ  Z 

save  CORRHA  CORRHA 
clear  CORRHA 
665 

'  loadDK; 
load  DM; 
load  DC; 

670  load  stat 

CORRHAD=[]; 

lempza=zeros(aset_size,aset_si2e); 

tempz=zeros(cset__size,cset_size); 

tempk=zeros(cset_size,cset_size); 

675  tempm=zeros(cset_size,cset_si2e); 

tempc=zeros(cset_size,cset_size); 
if  meters 

waitbarJiandle=waitbar(0,Tnstalling  DK,  DM,  &  DC...’); 
else 

680  dispC’Installing  DK,  DM,  DC) 

end 

for  i=  1  :uuinpoints-2 
if  meters 

waitbar(i/(iiumpoiiits-2)) 

685  end 

kcoiTected=kstat; 

mcorrected=mstat; 

ccorrected=cstat; 

690  tenip=DK(ndx3d([l  cset_size*(cset_size<-l)/2  (numpomts-(2))],l,l:cset_si2e*(cset  si2&fiy2  !))• 

skyindex=l;  ”  ’  ’ 

for  indexl=l;cset_size 

for  mdex2=indexl;cset_si2e 
tempk(indexl  ,index2)=temp(skyindex); 
tempk(index2,indexl  )=temp(skyindex); 
skyindex=skyindex+l ; 
end 
end 


00 


'05 


10 


5 


:o 


DK(:,1H]; 

temp“DM(ndx3  d([  1  cset_size  *(cset_size+ 1  )/2  (iiumpoiiits-(2 ))],!,!:  cset__size  *(cset  size+ 1  )/2  1 )  )* 
skyindex=l;  “  ’  ’ 

for  indexl=l:cset__size 

for  index2=indexi  :cset__si2e 
tempm(mdexl,mdex2)=temp(skyindex); 
teinpm(index2,index  1  )=terap(skyindex); 
skymdex=skyindex+l ; 
end 
end 

DMCUHJ; 

temp-DC(ndx3d([  1  cset_size*(cset_size-f- 1  )/2  (mimpoints-(2))],  1 , 1  :cset_size*(cset  size+I  )/2  1  ))■ 
skyindex=l;  “  /  5 

for  indexl=l:cset__size 

for  mdex2=indexl:cset__si2e 
tempc(indexl  ,index2)=temp(skymdex); 
tempc(index2,indexl  )=temp(skyindex); 
skyindex=skyindex+l ; 
end 
end 
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DC(;,1)=G; 

kcoiTected(cset_rel,cset_rel)=kcorrected(cset_reI,cset_rel)+tempk: 
mcorTected(cset_rel,cset_rel>=mcorrected(cset_rel,cset_rel)4-teinpin; 
ccorrected(cset_rel,cset_rel)=cx:orrected(cset  reLcset  relHlempc* 

725  ^  ^  F  . 

temp2a=kcorrected4-j*w(i+l)*ccoiTected-vv(i+l)^2*mcorTecled; 

tempha=mv(tempza); 

CORRIiAD=[CORRHAD  tempha(:)]; 
end 

730  if  meters 

close(vvaitbarJiandle) 

end 

clear  tempk  tempm  teinpc  lempza  tempha 
clear  kcorrected  mcorrected  ccorrected 
735  save  CORRHAD  CORRHAD 

clear  CORRHAD  DK  DM  DC 

%%%save  tlie 

Workspa<X%%%%%%%%%%%o/o%o/o%%%%%%o/o%%%%%o/o%%%o/oO/oO/aO/o%o/„o/^o/„o/^o/^o/^o/^o/^o/„o/^^^^^ 
740  save  INT 


745 
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clear 
closeall 
load  INT 
titles=0; 

5  ifpswitch==y 

wliitebg('wlute’) 

close 

end 

li=0; 

10  if  complete 

fignum=l; 
else 

fignum=21; 

end 

1 5  %%  write  frequencies  and  system  diag  to  diary  file  prt.out%%%%%%%%%%o^%o^o^o/^%o/^, 

ifexist('figOOO.ouf) 
delete  figOOO.out 
end 

diary  figOOO.out 
20  fprinlfC  \n’) 

fprintfC  \n') 
fprintf('A  setW) 
for  index=l:lengtli(aset) 

Q>rmtf([blaiks(4-Iengtli(mt2str(aset(mdex)))),mt2str<aset(mdex))l) 

25 '  ifrem(index,12)==0 
fprintf(V) 
end 
end 

fprintf(V) 

30  fprintf(V) 

fprintf('0  set\n') 

for  index=l  :lengtli(oset) 

fprintl([blanks(4-lengtli(int2str(oset(index)))),int2str(oset(index))]) 

if  reiu(index,  12)==0 
35  fprintf(V) 
end 
end 

tprintf(V) 

fprintf(V) 

40  fprintt'CComputed  C  set\n') 

for  index=l  :lengtli(cset) 

fpriiitf([blanks(4-leiigth(iut2slr(cset(mdex)))),int2s(i(cset(iiidex))]) 
if  rem(index,12)==0 
fprintf^V) 

45  end 
end 

fprmtf('\n') 
fprintf(V) 
fprintf(True  C  setW) 

50  for  index=  1  :lengtli(true_cset) 

fprintf([blanks(4-lengtli(int2str(tme_cset(index)))),int2str(true  cset(index))l) 
if  rem(index42)==0 
fprintfCVn’) 
end 

55  end 

fprintf(V) 

fprintf(V) 

bigger=max(lengtli(omegaa),lengtli(omegax)); 
if  complete 
60  allfreqs=[ ... 
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[omegaa/(2*pi);  zeros(bigger-lengtli(omegaa),l)] ... 

[omegax/(2*pi);  zeros(bigger-lengtIi(omegax),l)] ... 

else 

allfreqs=[ ... 

[oniegaa/(2*pi);  zeros(bigger-lengtli(omegaa),l )  ]... 

[omegax/(2*pi);  zeros(bigger-length(omegax),l)  ]... 
[omegared/(2*pi);zeros(bi^er-length(oinegared),l)  ] ... 

[omegaexto/(2*pi);  zeros(bigger-lengtli(oinegaexto),  1 )]. . . 

end 

[aaa  bbb]=size(allj5’eqs); 
fprmtf('System  Frequencies  (Hz)W) 
if  complete 

fprintfC  Anal  ExpW) 
bbb=2; 
else 

fprintfC  Anal  Exp  Red  Oset\n') 
bbb=4; 
end 

for  indexl=l:aaa 
prUine=[]; 
for  inde3{2=l:bbb 
if  allfreqs(indexl,index2)==0 
prtfreq=' 
else 

prtfreq=[blaiiks(8-lengtli(sprintf('%4g',allfreqs(indexl,mdex2)))),... 

spriiitf(%4g',allfreqs(indexl,index2))]; 

end 

prtline=[prtline,  prtfreq]; 
end 

prtline=[prtline,  Vi*]; 
fprintf(prtline) 
end 

diary  off 

%dos(’type  figOOO.out  >  Iptl  &'); 

xx=omegax(fmd(omegaa(lowinode)  <=  oinegax  &  omegax  <=  omegaa(liiglmiode)))/(2*DiV 
yx  1  =ones(length(xx),  1 ); 

aa=oniegaa(fuid(oinegaa(lowinode)  <=  omegaa  &  omegaa  <=  omegaa(lii2lmiode)))/(2*DiV 
ya  1  =ones(lengtli(aa),  1 ); 

red=omegared(find(omegaa(lo\vmode)  <=  omegared  &  omegared  <=  omegaa(higlunode))Vr2*r)iV 
yredl=ones(lengtli(red),  1 ); 

exto=omegaexto(find(omegaa(lowmode)  <=  omegaexto  &  omegaexto  <=  omegaa(lughniode))V('2*Di)' 
yextol  =ones(lengtli(exto),  1 ); 
if  complete 

xO=max([xx(lengtli(xx))  aa(lengtli(aa))]); 
else 

xO=max([xx(length(xx))  aa(length(aa))  red(length(red))  exto(lengtli(exto))]); 
end 

clear  omegax  omegaa  omegared  omegaexto 
clear  L  Z_ANAL  c_anal  c_exp  conn  li_exp  h_anal__red  k_anal  k_exp 
clear  kexta  kstat  kexto  m^anal  m_exp  niexto  mstat  temp J^diags 
if  complete 

mid__mdex=round(lengtli(cset)/2); 

else 

mid_index=round(lengtli(cset)/2); 

end 

load  H_ANAL_R 
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h=figure(h+l); 
figiiuni=fignum+l  ; 

plot(w/(2*pi),loglO(abs(H_ANAL_RED(ndx3d((l  aset_size»(asel_sizefl)/2  numpointsl ... 
1  ,sk>Ted(cset__rel(mid Jndex),cset_rel(niid_index)),V)))),’r~’,. 

125  w/(2*pi),loglO(abs(H_EXP(udx3d([l  aset__size*(aset_size+l)/2 

l,skyred(cset_rel(niidJndex),cset_rel(iiiidJndex)),V)))Xl3-.') 

v=axis; 

ylabel([  H(',mt2str(cset(niid_mdex)),V,mt2  str(cset(mid_mdex)),')  iii/lbf  (Log  1 0  of)'] ) 
xlabel(’0niega  (Hz)') 

130  if  titles 

title('Aiialytic  FRF  vs  Experimental  FRF‘) 
end 

hold  on 

ifabs(v(4)-v(3))<=  1 
135  v2(l)=v(l); 

v2(2)=v(2); 

v2(3)=v(3)-.  1  *abs(v(4)); 
v2(4)=v(4)+'.l*abs(v(4)); 
axis(v2); 

140  end 

grid  on 

lili=slegend(’ANALYnC','EXPERIMENTAL'); 

axes(lili); 

hold  off 

145  ifpsvvitch==y 
prtfig(fignum) 
delete(h) 
end 

150  clear  H^ANAL^RED 

clear  H_EXP 

if  complete 
h=figure(h+l); 

1 55  fignimi=figniun+ 1 ; 
else 

h=figure(h+2); 

fignum=fignum+2; 

end 

160  loadL 

plot(aset,abs(diag(L)),aset,abs(diag(L)),'*') 
ylabel('L(DOF)lbf/in’) 
if  complete 
xlabel(’DOF) 

165  else 

xlabel('ASETDOF) 
end 
grid  on 
if  titles 

1 70  title(['Localization  Matrix  Diagonal  at  Omega  =  ',nimi2str(\vfreq/(2*pi)), '  Hz'])- 

end  ’  ^  " 

if  pswitch==’y 
prtfig(figniun) 
delete(h) 

75  end 

clear  L 

h=figure(h+l); 
fignum=fignmn+ 1 ; 

BO  load  L^DIAGS 
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mesh(\v/(2*pi),aset,loglO(abs(L_DIAGS))) 

if  titles 

title('Frequeucy  Dependence  of  Localization  Matrix  Diagonals’) 
end 

1 85  xlabel(’Oniega  (Hz)’) 

ylabel(DOF) 

zlabe^T^DOF,  omega)  Ibf/in  (LoglO  oQ’) 
grid  on 

if  pswitch=y 
190  prtfig(flgnum) 

deiete(li) 
end 

h=figure(h+l); 

1 95  fignum=fignum+ 1 ; 

subplot(2,l,l) 

plot(w/(2*pi),log  I0(L_DIAGS(ndx3d([aset_size  1  nuinpoints],cset_rel(inid_index+2),  1 ,':')))) 
ylabel([’Error  Coord  ^;mt2s^^(cset(lludJndex+2));;,int2st^(cset(lnid__mdex4■2)),’)’]) 
if  titles 

200  title(’lbf7in  (LoglO  of)’ ); 

end 
grid  on 
vl=axis; 
subplot(2,l,2) 

205 

if  lengtli(cset_rel)  <  length(aset) 
lioldset=l:lengtli(aset) ; 
holdset(cset_rel)=[] ; 
index=round(lengtli(holdset)/2); 

210 

piot(w/(2*pi),loglO(L__DIAGS(ndx3d([aset_size  1  nunipoints],holdset(index),l ;;')))) 
ylabel(rNon-EiTor  Coord  L(’,iiit2str(aset(holdset(mdex)));;4nt2sti(aset(holdseKmdex)));)’]) 
if  titles 

title(Frequency  Dependence  of  L  Matrix  Non-Error  set  Diagonal  Elements' ) 

215  end 

v=axis; 
hold  on 
v2(l)=v(l); 
v2(2)=v(2); 

220  v2(3)=v(3>f((v(4)-v(3))».5Kvl(4)-vl(3)r.5 ; 

v2(4)=v(3>+<(v(4>v(3))*.5>f(vl(4)-vl(3))*.5; 
axis(v2); 
grid  on 
end 

225  xlabel('Omega  (PIz)') 
hold  off 
if  psvritch==’y' 
prtfig(fignum) 
delete(h) 

230  end 

clear  L_DIAGS 

loadZ__ANAL_R 
ioadZ_EXP 
235  h=figure(h+l); 

fignum=fignum+ 1 ; 

plol(w/(2*pi),ioglO(Z_ANAL_RED(ndx3d([l  aset_size*(aset_size+l  )/2 
nmnpoints],  1  ,skyred(cset_rel(mid_index),cset_rel(mid_index));:’)));r-’,. 
\v/(2*pi),logl0(Z_EXP(ndx3d([l  aset_size*(aset_size+l  )/2 
240  nmnpoints],  1  ,skyred(csel_rei(mid_iiidex),cset_rel(mid_index))/:'))),’b-.') 
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245 


250 


255 


260 


265 


270 


275 


280 
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!95 


00 


xlabel('0mega^2  (Hz)’) 

ylabel(pZ(',mt2str(cset(mid_mdex))/;,int2stit;cset(nudJiidex));^  IWin  (Loc^IO  of)’l) 

if  titles  ® 

title('Aiialytic  Impedance  vs  Experimental  Impedance') 
end 
v=axis; 
hold  on 

if  abs(v(4)-v(3))  <=  1 
v2(l)=v(l); 
v2(2)=v(2); 

v2(3)=v(3>.l*abs(v(4)); 
v2(4)=v(4>4-.I*abs(v(4)); 
axis(v2); 
end 
grid  on 

if  complete==  0 

hli=slegendCANALYnC  IMPEDANCE','EXPER1MENTAL  IMPEDANCE'); 
Ui=slegend('ANALYT[C  IMPEDANCE',EXPERIMENTAL  IMPEDANCE'); 

axes(lili)  - 

hold  off 
ifpswitch==y 
prtfig(fignmn) 
delet^h) 
end 

clear  Z^ANAL^RED  Z_EXP 


load  true_errs 
loadDK 
h=rigiire(h+l) ; 
fignuni=fignuni+l  ; 
if  complete 

^  plot(w(l:numpoints-2)/(2*pi),DK(ndx3d([l  cset_size»(cset_sizefl)/2 


numpoints-2],  1  ,skycsel(mid_  index, mid_index),';')),'r- 


Ml:numpoints-2)/(2*pi),tme_stiffiiess(cset(midJndex),cset(imd_index))*ones(l,nmnpoints-2)'b- 

ylabel([DK(',int2str(cset(mid_index)),',',mt2str(csel(mid  index)) ')  Ibf/in'l) 
xlabel('Omega  (Hz)’) 
else 

plot(\v'(l:numpoints-2y(2*pi),loglO(abs(DK(ndx3d([l  cset_size*(cset_size+-l)/2  numpoints- 
2],l,skycsel(mid_index,mid_index), ':')))), 'r-',... 


') 


w(l;nmnpoints-2)/(2*pi),logl0((true_stiffiiess(cset(mid_index),cset(mid_index))==0)+ 
true_stiffiiess(cset(niid_index),cset(mid_index))).*ones(l,nmnpoints-2),'b-.') 
ylabel([DK(',int2str(cset(mid_index)).',',int2str(cset(inid  index)),')  lbf/in’(Logl0  of)'l) 
xlabelCOmega  (Hz)')  "  y  o  ju 

end 


v=axis ; 
if  titles 

title('Computed  Stiffiiess  Error  vs  True  Stiffiiess  Error') 
end 

v2(I)=v(l); 
v2(2)=v(2); 
if  complete 
if  abs(v(4)-v(3))  <  1 
v2(3)=v(3 -hi  00*abs(v(4)-v(3)); 


v2(4)=v(4)+100*abs(v(4>v(3)); 

v(3)=v2(3); 

v(4)=v2(4); 

end 
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if abs(trae_stifrness(cset(midjndex),cset(midjndex))  -  v(3))  <  25*abs(v(4VvP)) 
v2(3)=v(3)-.5*abs(v(4>v(3)); 
v2(4)=v(4); 

elseif  abs(true_stifrness(cset(niid_mdex),cset(mid_index))  -  v(4))  <  .25*abs(v('4)-vn)) 

305  v2(4)=v(4)4*.5*abs(v(4>v(3)) ; 

v2(3)=v(3); 
else 

v2(3)=v(3); 

v2(4)=v(4); 

310  end 
else 

if  abs(loglO((true_stifriiess(cset(niid__index),cset(imd Judex))  ==  0)+... 
true__stiffiiess(cset(midjndex),cset(midjndex)))  -  v(3))  <  .25*abs(v(4>v(3)) 
v2(3)=v(3)-.5*abs(v(4>v(3)); 

315  v2(4)=v(4); 

elseif  abs(loglO((true__stiffiiess(cset(inidJndex),cset(mid Judex))  ==  0)+... 
true_stiffiiess(cset(inidjndex),cset(inid_mdex)))  -  v(4)7<  .25*abs(v(4)-v(3)) 
v2(4>=v(4>f.5*abs(v(4>v(3)) ; 
v2(3)=v(3); 

320  else 

v2(3)-v(3); 

v2(4)=v(4); 

end 

end 

325  axis(v2) 
grid  on 

Wi=slegendCCOMPUTED  STIFFNESS  ERROR’,TRUE  STjGFFNESS  ERROR’); 
axes(liii)  ; 

if  pswitch=='y’ 

330  prtfig(fignum) 

delete(h) 
end 

clear  DK 
load  DM 

335  h=figure(h+l)  ; 

fignuin=fignmn+l ; 
if  complete 

plot(\v(l:nuinpoints-2)/(2*pi),DM(ndx3d([l  cset_size*(cset_size+l)/2  numpoints-2],l,skycset(mid Judex, mid Jndex),’:’)),’r~ 

340  w(l:numpomts-2y(2*pi),true_mass(cset(nud_index),cset(mid_index)).*ones(l,nimipoints-2),'b-.') 
ylabel(['DM(',mt2str(cset(mid_index)),V,int2str(cset(inid_index)V)  lbf-sec''2/in']) 
xlabel('Oinega  (Hz)’) 
else 

pIot(w(l:numpomts-2)/(2*pi),loglO(abs(DM(ndx3d([l  cset^size*(cset_size+l)/2  numpoinls- 
345  2],  1  ,skycsel(midJndex,midJndex),V)))),’r-',... 

vv(  1  :numpomts-2)/(2*pi),log  1 0(... 

(-■• 

tnie_mass(cset(mid_index),cset(niid Judex))  ==  0  ... 
tme_mass(cset(mid Jndex),cset(mid Jndex)). . . 

350  ).  *ones(  1  ,numpoints-2),’b-.') 

yIabel([bM(',int2str(cset(mid_mdex)),V,iiit2str(cset(mid Judex)),’)  lbf-sec^2/in  (LoglO  of)'l) 
xlabel(’Omega  (Hz)’) 
end 
v=axis; 

355  if  titles 

title(’Computed  Mass  Error  vs  True  Mass  Error’) 
end 

hold  on 
v2(l)=v(l); 

360  v2(2)=v(2); 
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if  complete 
ifabs(v(4)-v(3))<  1 
v2(3)=v(3>10000*abs(v(4)-v(3)); 
v2(4)=v(4>fl0000*abs(v(4)-v(3)); 

365  v(3)=v2(3); 

v(4)=v2(4); 

end 

if abs(trae_niass(cset(mid_mdex),cset(niid_mdex))  -  v(3))  <  25*abs(vr4')-v('3')) 
v2(3)=v(3>.5*abs(v(4)-v(3)); 

370  v2(4)=v(4); 

elseifabs(true_raass(cset(mid_index),cset(mid_index))  -  v(4))  <  .25*abs(v(4W3')') 
v2(4)=v(4)+.5*abs(v(4>v(3)); 
v2(3)=v(3); 
else 

375  v2(3)=v(3); 

v2(4)=v(4); 

end 

else 

if  abs(loglO((true_mass(cset(mid_index),cset(mid_mdex))  ==  0)+... 

380  true_mass(cset(inid_mdex),cset(mid_index)))  -  v(3))  <  .25*abs(v(4Vv(3)) 

v2(3)=v(3)-.5*abs(v(4)-v(3)) ; 
v2(4)=v(4); 

elseif  abs(logl  0((trae_mass(cset(mid_index),cset(mid_mdex))  =  0>f... 
true_mass(cset(nud_index),cset(mid_index)))  -  v(4))  <  .25*abs(v(4Vvf3)) 

385  v2(4)=v(4)+.5*abs(v(4)-v(3)) ; 

v2(3)=v(3); 
else 

v2(4)=v(4); 

v2(3)=v(3); 

390  end 
end 

axis(v2) 
grid  on 

lili=slegend(’COMPUTED  MASS  ERROR’,TRUE  MASS  ERROR’)- 

395  a\es(lili) 

if  pswitch— y 
prtfig(fignum) 
delet€<h) 
end 

400  clear  DM 

load  DC 
h=figure(h+l); 
fignuni=fignuin+ 1 ; 

105  if  complete 

plot(w(l:numpomts-2y(2*pi),imag(DC(ndx3d([l  cset_size*(cset_size+I)/2  numpoints- 
2],l,skycset(mid_index^iid_index),';'))),'r-',... 

\^l:nmiipomts-2y(2»pi),iinag(true_dainpiiig(cset(mid_index),cset(imd_mdex)).*ones(l4iumpoints-2))'b-') 
ylabel(['DC(',mt2str(cset(nud_index)),V,int2str(cset(mid_index)),')  M-sec/mH) 

10  xlabel('Oniega  (Hz)') 
else 

plot(w(l:numpomts-2y(2*pi),loglO(abs(DC(ndx3d([lcset_si2e*(cset_size+iy2numpoints- 

2],l,skycset(mid_index^iid_mdex).’:')))),'r-',... 
w(l  :numpoints-2y(2*piXloglO(ab<... 

15  (... 

true_damping(cset(mid_index),cset(mid_index))  =  0 ... 

)+true_daniping(cset(mid_index),cset(mid_mdex)). . . 

).  *ones(  1 4iumpoints-2)),'b-.') 

ylabel(['DC(',mt2str(cset(mid_index)),V,int2str(cset(mid  index)),’)  lbf-sec/in(Logl0of)'l) 

20  xlabelCOmega  (Hz)')  “  ,  v  e  /j; 
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end 
v=axis; 
if  titles 

title('Coinputed  Damping  Error  vs  True  Damping  Error’) 

425  end 

v2(l)=v(l); 
v2(2)=v(2); 
if  complete 
ifabs(v(4)-v(3))<  1 

430  v2(3)=v(3)-100*abs(v(4)-v(3)); 

v2(4)=v(4>fl00*abs(v(4).v(3)); 

v(3)=v2(3); 

v(4)=v2(4); 

end 

435  if  abs(true_damping(cset(midjudex),cset(mid_index))  -  v(3))  <  .25*abs(v(4)-v(3)) 
v2(3)=v(3>.25*abs(v(4>v(3)) ; 
v2(4)=v(4); 

elseif  abs(true_dampmg(cset(mid_index),cset(midjndex))  -  v(4))  <  .25*abs(v(4)-v(3)) 
v2(4)=v(4)+.25*abs(v(4)-v(3)) ; 

440  v2(3)=v(3); 

else 

v2(3)=v(3); 

v2(4)=v(4); 

end 

445  else 

if  abs(loglO((true_damping(cset(mid_mdex),cset(midJndex))  ==  0>i-... 
true_dainping(cset(iTud_^mdex),cset(nnd__index)))  -  v(3))  <  .l*abs(v(4)-v(3)) 
v2(3>=(v(3>.25*abs(v(4)-v(3))) ; 
v2(4)=v(4); 

450  elseif  abs(log  1 0((true__dampiiig(cset(mid_index),csel(mid_index))  ==  0>f. . . 

true_damping(cset(mid__index),cset(mid_index)))  -  v(4))  <  .l*abs(v(4)-v(3)) 
v2(4)=(v(4)+.25*abs(v(4)-v(3))) ; 
v2(3)=v(3); 
else 

455  v2(3)=v(3); 

v2(4)=v(4); 

end 

end 

axis(v2); 

460  grid  on 

Wi=slegend('COMPUTED  DAMPING  ERROR’,’TRUE  DAMPING  ERROR’); 
axes(liii); 
ifpswitch=='y’ 
prtflg(fignum) 

465  delete(h) 
end 

clear  DC  true_damping  true__mass  true_stiffiiess 

470  %%%% 

loadH__EXP 

loadCORRHA 

%for  mid__iudex=  1  :min(cset__si2e,multplot) 

475  h==rigure(h+l); 

figiuun=fignum+ 1 ; 

piot(\v/(2  *pi),log  1 0(abs(H_EXP(ndx3d(  [  1  aset_size*(aset_size+ 1  )/2 
numpoints],  1  ,skyred(cset_rel(niid>idex),cset_rel(mid__mdex));:’))));r--’,. . . 

\v/(2*pi),loglO(abs(CORRHA(ndx3d([aset_size  aset^size  numpoints],cset_rel(mid_index),cset_rel(midJndex);:’)))),V.’) 

480  ylabel(fH(’,mt2str(cset(mid^index));;,int2str(cset(mid_index));)  in/lbf  (LoglO  of)']) 
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xlabel(’Oiiiega  (Hz)’) 
if  titles 

titleCExperiinental  FRF  vs  DZ  Corrected  FRF') 
end 

485  grid  on 

v=axis; 

Wi=slegend(’EXPERIMENTAL  FRF/CORRECTED  FRF); 
axes(lili)  ; 
if  pswitch=='y’ 

490  prtfig(fignuni) 

delete(h) 
end 

clear  CORRHA 
495  load  CORRHAD 

h=ligure(h+l); 
flgniim=fignuin+l ; 

plot(  w/(2  *pi),log  1 0(abs(H_EXP(iidx3d([  1  aset_si2e*(aset_size+ 1  )/2 

nunipoints],l,skyred(cset__rel(midJiidex),cset__rel(midJndex));:'))));r--’,... 

500  w(  1  :niunpoints-2)/(2*pi),Iog  1 0(abs(CORRPIAD(ndx3d([aset_size  aset_size  (numpoints- 

2)],cset_rel(mid_index),cset__rel(niid_index),':')))),’b-.’) 
ylabel([’H(',int2str(cset(midJndex));;,int2str(cset(midJndex));)  in/lbf  (LoglO  oO'l) 
xlabel(’Qinega  (Hz)’) 
if  titles 

505  title(’Experiinental  FRF  vs  DK/DM/DC  Corrected  FRF’) 

end 
grid  on 
v=axis; 

lili=slegend('EXPERIMENTAL  FRF',’CORRECTED  ANALYTIC  FRF’); 

510  axes(lili); 

if  pswitch=’y' 
prtfig(fignmn) 
delete(h) 
end 

515  clear  CORRHAD 

save  plotniun  fignuni 
clear 


103 


CHAP3.M 


%%%%SSTCONF  FILE  FOR  A  SPATIALLY  COMPLETE  BEAM 
beaiiinidl  ;%beamnidl  is  a  cantilevered  20  element  beam  model 
%siimndl  ;%sinimdi  is  a  masses  and  springs  model 
pswitch=’n’  ;%%do  we  print? 

5  titles=l  ;%%display  titles 

meters=l  ;%use  progress  meters 
whitebg(l3lack’)  ;%switch  to  black  figure  background 
close 

lowmode=l; 

10  luglmiode=10; 

static=l;%%%%%%%%%%%%%reduction  metliod  0=Guyan,l=IRS,2=Extraction 
%%define  A  set  &  O  set%%%%%%%%%%%%%%%%%%%o/oO/oO/^o/^o/^,o/^o/^^o/^o^o/^^^^ 
aset=l:numdof; 

1 5  %oset=2:2:numdof;  %%A  SPATIALLY  INCOMPLETE  BEAM%%%%%%%%%% 

%OSet=SOrt(oset);%%%%%%%%%%%%%%%%o/o%%o/oO/oO/oO/o0^o/^o/^o/^o/^o/^o/^^^o/^o^^^^^ 

oset=[];  %%%%A  SPATIALLY  COMPLETE  BEAM%%%%%%%%%%%%o/o%o/o 
20  save  sstconf 
sst 

plotsst 
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10 


%%%%SSTCONF  FILE  FOR  A  SPATIALLY  INCOMPLETE  BEAM 

beanundl  ;%beaiiiiiidl  is  a  cantilevered  20  element  beam  model 

%simmdl  ;%simmdl  is  a  masses  and  springs  model 

pswitch=’n’  ;%%do  we  print? 

titles=l  ;%%display  titles 

nreters=l  ;%use  progress  meters 

wlutebg(’black’)  ;%svvitch  to  black  figure  background 

close 

lowmode=l; 

liiglimode=10; 


slatic=l;%%%%%%%%%%%%%reduction  metiiod  0=Guyan,l=IRS,2=Extraction 
%%defme  A  set  &  O  set%%%%%%%%%%%%%%%%%%%%%%o/o%%%%%%%%%o/o% 
aset=l:numdof; 

1 5  oset=2:2:nuiudof;  %%A  SPATIALLY  INCOMPLETE  BEAM%%%%%%%%%% 

oset=sort(oset);  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%o^% 


save  sstconf 
sst 

20  plotsst 


beammdl  ;%beanundl  is  a  cantilevered  20  element  beam  model 
%simmdl  ;%simmdl  is  a  masses  and  springs  model 
pswitch=’n’  ;%%do  we  print? 

25  titles=l  ;%%display  titles 

meters=l  ;%use  progress  meters 
wMtebg(l3lack’)  ;%switch  to  black  figure  background 
close 

lowmode=l; 

30  liiglunode=10; 


%%define  A  set  &  O  set%%%%%%%%%%%%%%%%%%%%o/o%%%%o/o0/o0/^o/j,o^o^o^o/^o^ 
aset=l:numdof; 

oset=2:2:nuindof;  %%A  SPATIALLY  INCOMPLETE  BEAM%%%%%%%%%% 

35  oset=sort(oset);%%%%%%%%%%%%%%%%%%%%%%%%%o/o%o/oO/oO/oO/o0^o^o^o/^%o/Qo^o/^^ 


static=2; 
save  sstconf 
dofig4_3 
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%%%%%%%%%%%%%%%®^%0/O0^0^%%0^0^0/J0/„%0^0^0^%0/^,0/„0^0^%0^0^%0^0^0^0/„0^0^0^0^0^%0^0^0/^0^%^^^ 

%%%%%%%%%%%o/o 

%%%  COMPARE  SINGLE  MODE  MATRIX  AND  INTEGRAL  %%%%%%%%%%%%%%«^®/o%%%%%%%%%%%% 
%%%  SOLUTIONS  UNDER  VARIOUS  CONDITIONS  %%%%%%%%%%%%%%%%%%%%<>/o%%%%o/o%o/g 

6i%%%%%%%%%%% 

clear 

closeall 

loadINT 

1 0  whitebg(\vhite’) 

clear  L  Z__ANAL  c__anal  c_exp  comi  h_exp  li_aiial_red  k_aiial  k_exp 
clear  kexta  kexto  m__aiial  m_exp  mexta  mextO  tenip_l_diags 
clear  true^dainping  tmejnass  true^stiflhess  z_anal  z_anal__red  z_exp 
clear  L  Z_ANAL  c_anal  c_exp  comi  h_exp  h_anal_red  k_anal 
1 5  clear  kexta  kexto  m_aiial  mexta  mexlO  temp_l__diags 

loadH__EXP 
load  exp 

cset=[5  7  9  11  13  15] 

20  cset_^rel=[34  5  6  78] 

cset_size=leiigtli(cset) 
mid_mdex=round(lengtli(cset)/2); 
fuieiaess=200 
ODZl=[]; 

25  tempi  =[]; 

odfreq=omegax(l); 

center_freq=odfreq/(2*pi) 

wl=odfreq; 

odiength=2*pi; 

30  odivisioiis=2; 

lowerfreq=odfreq- 1 5  *  odi  engtli 
upperfreq=odfreq+25  *odlengtli 

ow  1  =odfi'eq-.5  *odlengtli:odleiigtli/odi  visions:odfreq+.  5  *odleiigtli; 
partition=o\vl/(2*pi) 

35  df=(odlengtli/odivisions)/(2*pi) 

for  couiit=  1  ;lengtii(ovvl ) 
z_anal_red=kstat-owl(couiit)^2*mstat; 
h_aual_red=inv(z_anal_red); 

z_exp=k_exi>fj  *owl  (couiit)*c_exp-ovvl  (count)'^2  ♦m^exp; 

40  h_exp=inv(z__exp); 

h_exp=h_exp(aset,aset); 

hacc=h_aiial_red(cset_rel,cset_rel); 

lixcc=li_exp(cset_rel,cset_rel); 

dz=inv(iiiv(mv(hacc)*(hacc-lixcc)*mv(hacc))-hacc); 

45  ODZl=[ODZl;dz]; 

templ=[templ;  eye(cset_size)  -owl(coimt)^2*eye(cset_size)  +j*owl(coimt)*eye(csel_size)]; 
end 

%aieleftdz=ODZl 
50  %tliebigone=templ 

DKDMDC 1  =templ\ODZl ; 

ODK 1  =DKDMDC  1  ( 1  :cset_size, : ) 

ODMl  =DKDMDC  1  (cset_size+ 1 :2  *cset_size,;) 

ODC 1  =DKDMDC  1  (2*cset_size+ 1 :3  *cset_size,:) 

55  save  ODMl  ODMl  DKDMDCl  ODZl  owl  wl 

saveODKl  ODKl 
save  ODC  1  ODCl 

%-HH  H  I  1  I  I  I  1  1  I  I  1  i  I  1  I  I  l  i  11  i  I  [  I  I  1  I  j  I  I  I  I  M  WWW  i  }  1  I  I  M  I  M  I  1  1  I  1  [  I  i  I  I 

60  load  stat; 
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odcorrlia=[]; 

w=lowerfreq:(upperfreq-lowerfreq)/fiiieiiess:upperfreq; 

nuinpoints=lengtli(w); 

kcorrected=kstat; 

65  mcorrected=nistat; 

ccoiTected=cstat; 

kcorrected(cset_rel,cset_rel)=kcorrected(cset_rel,cset_rel>fODK  1 ; 
mcorrected(cset_rel,cset_rel)=incon‘ected(cset_rel,cset_rel)+ODM  1 ; 
ccorrected(cset_rel,cset_rel)=ccoiTected(cset_rel,cset_rel>^^I^  1 ; 

70  for  i=l  aiuinpoiiits 

teinpza=kcoiTected+j  *w(i)*ccoiTected-w(i)^2  *mcoiTected; 
tempha=mv(tempza); 
odcorrha=[odcorTha  tenipha(:)]; 
end 
75 

uncorrha=[]; 
for  i=l:nunipoints 

tempza=kstat+j*w(i)*cstat-w(i)^2*instat; 

tempha=inv(tempza); 

80  uncorrha=[uncorTha  templia(:)]; 
end 

if  meters 

waitbar_handle=waitbar(0;Computing  Experimental  FRF'); 

85  else 

dispCGetting  experimental  FRF') 
end 

H_EXP=[]; 

for  count=  1  :length(w) 

90  if  meters 

waitbar(count/lengtli(w)); 

end 

2_exp=k_expfj  *  w(coimt)*c_exp-\v(count)^2*m_exp; 
h_exp=mv(2_exp); 

95  li_exp=h_exp(aset,aset); 

skymdex=l; 
for  indexl = 1  :aset__size 
for  index2=indexl  :aset_si2e 
redJiolder(skymdex)=h__exp(indexl  ,index2); 

1 00  skyindex=skyindex+l ; 

end 
end 

H_EXP=[H_EXP  red^holder’]; 
end 

1 05  if  meters 

cIose(waitbar_handIe) 

end 

clear  templia  tempza  kcorrected  mcorrected  ccorrected 
110  figure(l); 

plot(w/(2  *pi),log  1 0(abs(H_EXP(ndx3d([  1  aset_si2e*(aset_si2e+'  1  )/2  iiumpoints],. 

1  ,skyred(cset__rel(inid_index),cset_rel(mid_index));:'))));r-',. , . 

w/(2*pi),logl  0(abs(odcorrha(ndx3d([aset_size  aset_size  nuinpoints],cset_rel(inid__index),cset_rel(niid Jndex);:'))));g--V 

w/(2*pi),loglO(abs(uncoiTha(ndx3d([aset_sizeaset_sizenuinpomts],cset_rel(mid_mdex),cset_rel(mid_mdex),':')))),'b:') 

ylabel([’H(',mt2str(cset(mid_mdex)),V,mt2str(cset(imd_mdex)),')  iii/lbf(logl0  oO']) 

20  xlabeI('Omega  (Hz)') 
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%title('Single  Mode  Corrected  FRF  (Mode  1 ) ') 
grid 

lili=slegend('Experimental  FRF,  Mode  1  Corrected  FRF',  'Uncoirecled  Analytic  FRF'); 

125  axes(lili) 

print  -dnifile  fig5_l 


130 


135 


140 


145 


into  DM,  DK,  &  DC  %%%%%%%%%%%%%%%%%%%%%%o^o^%%^^^^^y^^^y^y^ 
%%%  USING  MATEUX  TECHNIQUE  %®^%%%%%%%%%%%%%%%%%%%%%%%%o^o^%%%o^o^o/^JO/^J 

%%%%%%%%%%%«/oO/oO/oO/„0/„0/,0/„0/,0/„0/„0/„0/„»/„0/„0/„0/„0/„0/„0/„0/„0/„0/„0/„0/„«/,0/„0/„0^^^^^ 

WGRi^J^%%%%%%% 


templ=[]; 

odlengtli=2*pi; 


owl-omegax(l)-5♦odlengtll;(omegax(2)-omegax(I>^-10*odIengal)/14:omegax(2>f5♦odlellgtll; 


for  counl=  1  :lengtli(owl ) 
z_anal_red=kstat-owl(count)^2*instat; 
h__anal_red=inv(z_anal_red); 

z_^exp=k_exp+j  ♦owl  (count)*c_exp-owl  (couiit)^2  ♦m_exp; 

h_exp=inv(z_exp); 

h_exp=h__exp(aset,aset); 

hacc=h_anal_red(cset_rel,cset_rel); 

lixcc=h_exp(cset_rel,cset_rel); 

dz=mv(inv(mv(hacc)^(liacc4ixcc)*mv(hacc)>.hacc)' 

ODZl=[ODZl;  dz]; 

templ=[templ;  ey<cset__size)  -owl(count)^2*eye(cset_size)  +j*owl(count)*eye(cset_^size)]; 


end 

150 

DKDMDCl=templ\ODZl; 
ODKl=DKDMDCl(l:cset_size,:) ; 

ODMl  =DKDMDC  1  (cset_size-f“  1 :2  ♦cset_size,:) ; 
ODC 1  =DKDMDC  1  (2*cset_size+ 1 :3  ♦cset_size,: ); 
155  save  ODMl  ODMl  DKDMDCl  ODZl  owl  wl 
saveODKlODKl 
saveODCl  ODCl 


160 


165 


170 


175 


180 


odcorrlia=[]; 

w=omegax(l>15*pi:(omegax(2)-omegax(l>f35*pi)/rmeness:oinegax(2)+20*pi; 

nuinpoints=lengtli(w); 

kcorrected=kstat; 

mcorrected=mstat; 

ccorrected=cstat; 

kcorrected(cset_rel,cset_rel)=kcorrected(cset_rel,cset_rel>K)DKl; 
mcorrected(cset_rel,cset_rel)=mcorrected(cs^__rel,cset__reI>+ODMl; 
ccorrected(cset_rel,cset_rel)=(xorrected(cset_rel,cset_rel>^DC  1 ; 

for  i=l;nuinpoints 

tempza=kcoiTected+j*w(i)*ccorrected-w(i)^2^mcorrected; 
tempha=inv(tempza); 
odcorrha=[odcorrha  tempha(:)]; 
end 

uncoiTlia=[]; 
for  i=l:niimpoints 

tempza=kstat4j^w(i)*cstat-w(i)^2*nistat; 
tempha=inv(tempza); 
uncorrha=[uncon‘ha  tempha(;)]; 
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end  ' 
if  meters 

wailbarJiandle=waitbar(0;Coniputing  Experimental  FRF')- 
185  else 

disp('Getting  experimental  FRF) 
end 

H_EXP=0; 

for  comit=  1  :lengtli(w) 

190  if  meters 

waitbar(count/leiigtli(w)); 

end 

z_exp=k_exp+j  *vv(comit)*c__exp-w(count)''2  *m_exp; 
li_exp=inv(z_exp); 

1 95  h_exp=h_exp(aset,aset); 

skyindex=l; 
for  index  1 = 1  :aset_size 
for  mdex2=indexl;aset_size 
red_holder(skyindex)=h_exp(indexl  ,index2); 

200  skyindex=skyindex+l ; 
end 
end 

H_EXP=[H_.EXP  redjioldei^; 
end 

205  if  meters 

close(waitbarjiaiidle) 

end 

clear  tempha  teinpza  kcorrected  mcorrected  ccorrected 
210  rigure(2); 

plot(w/(2*pi)doglO(abs(H_EXP(ndx3d([l  aset_size*(aset_size+l  )/2  numpoiiits],... 
l,skyred(cset__rel(mid Judex), cset_rel(mid__index));:'))));r-',... 

w/(2  *pi),log  1 0(abs(odcorrha(ndx3d([aset_size  aset_size  numpoints],cset__rel(mid Jndex),cset_rel(niid  index),':’)))),'g~',. 

w/(2*pi),loglO(abs(uncoiTha(ndx3d([aset_sizeaset_sizenumpoinls],cset_rel(inidJndex),cset_rel(mid_index),';')))),*b:’) 
ylabel(['H(’,int2str(cset(midJndex)),';,int2str(cset(midJndex)),')  in/lbf (loglO  of)']) 

220  xlabelCOmega  (Hz)') 

%title('Single  Mode  Corrected  FRF  (Mode  1) ') 
grid 

Wi=slegend(Experimental  FRF',  'Single  Mode  Corrected  FRF,  'Uncorrected  Analytic  FRF); 

125  axes(lili) 

print  -dmfile  tig5_2 


^30 


odivisions=3; 

for  ind=l  :lengtli(lenctrl) 

35  odlengtli=lenctrl(ind)*2*pi; 

ODZ=[]; 

teinp=[]; 

o\v=odfreq-.5*odlengtlnodlengtli/odivisions:odfreq+.5*odlength; 

for  counl=  1  :lengtli(ow) 

=10  z_ana  Jred=kstat-ow(comU)^2  *mstat; 


%%%%%%%%%0/oO/oO/oO/,0/,0/,0/,0/^0/,0/,0/,0/,0/,^^^^^^^ 

into  ODM,  ODK,  &  ODC  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%«/o%o/„o/o%o/oO/o«/oO/o./,o/„o/„o/,./,o/.o/,o/,o/„«/,o/„o^^^^^^^ 

pc&f^(^r%M^4fi^^(^1^fi&eq-lowerffeq)/fineness:upperfreq; 
lenclrl=[25  10  I  j; 
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h_aiial_red=iiiv(z_anal_red); 

z_exp=k_exp+j*ow(couiit)*c_exp-ow(coiuit)''2*m_exp; 

h_exp=inv(z_exp); 

h_exp=h_exp(aset,aset); 

245  hacc=h_anal_red(cset_rel,cset_rel); 
hxcc=h_exp(cset_rel,cset_rel); 
dz=inv(inv(inv(liacc)*(hacc-hxcc)*mv(hacc))-liacc); 

ODZ=[ODZ;  dz]; 

temp=[temp;  eye(cset_size)  -ow(count)''2*eye(cset_size)  +j*ow(count)*eye(cset_size)]; 
end 

DKDMDC=temp\ODZ; 

ODK=DKDME)C(l:cset_size,:) ; 

255  0DM=DKDME)C(cset_si2ei- 1 :2  *cset_size,:) ; 

0E)C=DKDMDC(2  *cset_size+ 1 :3  *cset_size,: ); 
clear  hacc  lixcc  Ii_exp  z_exp  z__anal_red  h_aiial__red  dz 
save  ODM  ODM  DKDMDC  ODZ  ow 
saveODKODK 
260  save  ODC  ODC 

waitbar__handle=waitbar(0;Computing  Experimental  FRF); 
else 

265  dispCGetting  experimental  FRP) 
end 

H_EXP=a; 

for  comit=l:lengtli(plotw) 
if  meters 

270  waitbar(count/lengtli(plotw)); 

end 

z_exp=k_exp+j*plotw(count)*c_exp-plotw(count)'^2*m__exp; 

li_exp=inv(z__exp); 

h_exp=li_exp(aset,aset); 

275  skyindex=l; 

for  index  1 = 1  :aset_size 
for  index2=indexl:aset_si2e 
red Jiolder(skyindex)=h_exp{index  1  ,index2); 
skyindex=skyindex+l ; 

280  end 

end 

H^EXP=[H_EXP  red^holdef]; 
end 

if  meters 

285  close(\vaitbar_handle) 

end 

odcorrha=[]; 

kcorrected=kstat; 

290  mcorrected=mstat; 

ccoiTected=cstat; 

kcorrected(cset__rel,cset_rel)=kcorrected(cset_rel,cset__rel>K)DK; 

mcorrected(cset_rel,cset_rel)=mcon-ecte^cset_rel,cset_rel>K)DM; 

ccorrected(cset_rel,cset_rel)=ccorrected(cset_rel,cset_rel>K)DC; 

295  for  i=  1  :Iengtli(plotw) 

tempza=kcorrected+j*plot\v(i)*ccorTected-piotw(i)'^2*mcoiTected; 
tempha=inv(tempza); 
odcorrha=[odcorrIia  tempha(:)]; 
end 

300  ifind==l 
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pIotl=odcoiTha; 

el$eifmd==2 

plot2=odcorrha; 

elseifmd==3 

plot3=odcoiTha; 

else 

plot4=odcorrha; 

end 

end 

figure(3); 

plot(plolw/(2*pi),loglO(abs(H_EXP(ii(i\3d([l  aset_size*(aset_size+l)/2  lengdi(plot\v)],... 
l,skyred(cset_rel(mid_mdex),cset_rel(raid_index)),';')))),'r-’,... 

^  plotw/(2  pi),  Iogl0(abs(plotl(ndx3d([aset_size  asel_size  lengtIi(plotvv)],cset_rel(mid_index),cset_rel(mid_index),':')))),'g— 
»••• 

plotw/(2*pi),  Iogl0(abs(plot2(iidx3d((aset_size  aset_size 

lengtli(plotw)},cset_rel(imd_mdex),cset_rel(mid_uidex),';')))),'b:',... 
plotw/(2*pi),  Iogl0(abs(plot3(ndx3d([aset_sizeaset_sizelengtli(plotw)],cset_rel(iiiid_mdex),cset_rel(mid_iiidex),';')))),'k-.') 

ylabel([H(',mt2str(cset(imd_index)),V,mt2str(cset(mid_mdex)).')  iii/lbf  (log  1 0  of)']) 

xlabel('Omega  (Hz)') 
grid  on 

"/otitleCEXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  5  SOLUTIONS') 
lili=legend('Exp  FRF, ... 

[nimi2str(lenctrl(l)),'  Hz  Bandwiddi  Corrected  FRF'], ... 

[iiuiii2str(lenctrl(2)),'  Hz  Bandwiddi  Corrected  FRF'], ... 

[nuni2str(lenctrl(3)),'  Hz  Bandwiddi  Corrected  FRF]) ; 
axes(hli) 

print  -dnifile  fig5_3 

odlengtii=l*2*pi; 
countctrl=[3  10  50  200] ; 
for  md=l:lengtii(countctrI) 
odivisions=countctrl(md)-l ; 

ODZ=[]; 

teinp=[]; 

ow=odfreq-.  5*odlengtli:odlengai/odivisions:odfreq+.5  *odIeiigtli; 
for  count=  1  :lengtli(ow) 

2_anal_red=kstat-ow(count)^2*nistat; 

h_aiial_red=mv(z_anal_red); 

z_exp=k_exp+j  *ow(count)*c__exp-o  w(count)'"2  *m_exp; 

h_exp=mv(z_exp); 

h_exp=h__exp(aset,aset); 

hacc=h_anal_red(cset_rel,cset_rel); 

Iixcx=h_exp(cset_rel,cset_rel); 
dz=mv(mv(inv(hacc)*(lia(x4L\(x)*inv(hacc)>hacc)- 
ODZ=[ODZ;  dz]; 

temp— [temp,  eye(cset_si2e)  -ow(count)^2*eye(cset__size)  +j*ow(count)*eye(cset  size)]; 
end 

DKDMDC=temp\ODZ; 

ODK=DKDMDC(l:cset_size,:) ; 

ODM=DKDME)C(cset_size+l  :2*cset__size,:) ; 

ODC=DKDMDC(2  *cset_size+ 1 :3  *cset_size,; ); 
clear  hacc  lixcc  h_exp  z^exp  z_anal_red  h_aiial_red  dz 
save  ODM  ODM  DKDMDC  ODZ  ow 
save  ODK  ODK 
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365 


370 


375 


380 


385 


390 


395 


405 


410 


415 


420 


save  ODC  ODC 

odcorrha=[]; 

kcoiTected=kstat; 

incx)jnrected=mstat; 

ccorrected=cstat; 

kcorrected(cset__rel,cset_rel)=kcon'ected(cset_rel,cset_reI>fODK; 

mcorrected(cset_rel,cset_rel)=mcoiTected(cset_rel,cset_rel)+ODM: 

cco^ected(cset_rel,cset_rel)=ccorrected(cset_rel,cset_rel>K)DC; 

for  i=  1  :length(plotw) 

teiiipza=kcorrected+j  *plolw(i)*cconrected-plot\v(i  )^2  *mcoiTecled; 
lempha=mv(tempza); 
odcoiTha=[odcorTha  tempha(:)]; 
end 

if  iud==l 
plotl=odcoiTlia; 
elseif  ind==2 
plot2=odcorriia; 
elseif  ind==3 
pIot3=odcorTha; 
else 

plot4=odcorrlia; 

end 

end 

figure(4); 

plol(plotw/(2*pi),logl0(abs(H_EXP(ndx3d([I  aset_size*(aset_size+l)/2  length(plotw)].... 

l,skyred(cset_rel(midJndex),cset_rel(niidJndex)),V)))),'r-',.-. 

^  Plotw/(2*pi),loglO(abs(plotl(ndx3d([aset_sizeaset_sizelengUi(plotw')],cset_iel(mid_index),cset_rel(mid_mdex),':’)))),'g- 

plotw/(2*pi),  Iogl0(abs(plot2(ndx3d([aset_size  aset_size 

lengtli(plot\v)],cset_rel(mid_index),cset_rel(niid_mdex),':')))),'b:',... 

^  Plotw/(2*pi)Jogl()(abs(plot3(ndx3d([aset_sizeaset_sizeleiigth(plotw)],cset_rel(midJndex),cset_rel(niidJndex),':')))),’k- 

. 

plot\v/(2*pi),logl0(abs(plot4(ndx3d([aset_sizeasel_sizelength(plolvv)],cset_rel(niid_iudex),cset_rel(mid_index),':'))));m:’) 

ylabel(['HC,int2str(cset(mid_index)XV,int2str(cset(inid_index)),')  iii/lbf(loglO  oOl) 

xlabel('Omega  (Hz)') 
grid  on 

’/otitleCEXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  I  SOLUTIONS’) 

Iili=legend('ExpFRF.... 

(nuin2str(countctrI(l)),'  Rs  Corrected  FRF'], ... 

[nuiH2str(couiitctrl(2)),'  Pts  Corrected  FRF], ... 

(iimn2str(couiitctrl(3)),'  Pts  Corrected  FRF'], ... 

[num2str(countctrl(4)),'  Pts  Corrected  FRF]); 
axes(lili) 

print  -dinfile  fig5_4 


%plotvv=lovv^erfreq:dAv:upperfreq; 
lenctri=[25  10  1]; 
odivisions=3; 
for  ind=l:lengtli(lenctrl) 
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odlengai=lenctrl(iiid)*2*pi; 

owl=odfreq-.5*odlengai;odlengtli/(odivisions):odfreq+.5*odlengtli; 

Wl=owl; 

%%compute  integrals  %%%%%%%%%%%%%%%%%%o/oO/o%o/„o/oO/oO/o%%o/„o/„o/^^^^^^ 
425  %weiglit=ones(si2e(owl  ))*(ome®ax(  1 )); 

weiglit=owl; 

DZ_R=[]; 

DZJ=[]; 

for  couiit=l  :IeiigtIi(owl ) 

430  z_anal__red=kstat/(j  *owl  (count))+j  *owl  (count)*mstat; 

h_anal_red=inv(z_anal_red); 

Z__exp=k^exp/0*owl(comit)>fc__exp/(j»owl(comit))+j*ovvl(couiit)*m  exp; 

h_exp=iiiv(z_exp); 

h_exp=h_exp(aset,aset); 

435  hacc=h_anal_red(cset_rel,cset_rel); 

lLxcc=h_exp(cset_rel,cset_rel); 
dz=mv(mv(mv(hacc)*(liacc-hx(x)*inv(hacc)H^^ 
dz_r=real(dz); 
cizj=miag(dz); 

440  DZ_R=pZ^Rdz^i(:)]; 

DZJ=[DZ J  dz_i(:)]; 
end 

[imXINTMINTC]=intsub(DZ_I,DZ  R,owl,Wl,l./(j*ow]),cset  size)- 
445  odcorrha=(];  "  v 

load  stat 

kcorTected=kstat; 

mcorrected=nistat; 

ccoiTected=cstat; 

450  kcoiTected(cset_rel,cset_rel)=kcorrected(cset_rel,cset_rel>fINTK; 
mcoiTected(cset_rel,cset_rel)=mcorrected(cset_rel,cset_rel>fINTM; 
cco^ected(cset__rel,cset_rel)=ccorrected(cset_rel,cset_rel)+l]SrrC; 

for  i=  1  :lengtli(plotw) 

tempza=kcorrected+j*plotw(i)*ccoiTected-plotw(i)^2*mcoiTected; 

455  tempha=inv(tempza); 

odcoiTha=[odcorrha  tempha(:)]; 
end 

ifmd==l 

intplot  1  =odcorrha; 

^0  savelNTMl  INTM 

save  INTKl  INTK 
save  INTCl  INTC 
elseif  ind==2 
mtplot2=odcoiTha; 

^65  save  INTM2  INTM 

save  INTK2  INTK 
save  INTC2  INTC 
elseif  ind==3 
intplot3=odcorrha; 

70  save  INTM3  INTM 

save  INTK3  INTK 
save  INTC3  INTC 
else 

intplot4=odcoiTlia; 

75  save  INTM4  INTM 

save  INTK4  INTK 
save  1NTC4  INTC 
end 
end 
50 
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if  meters 

waitbar_handIe=waitbar(0,'Computmg  Experimental  FRF); 
else 

disp(’Getting  experimental  FRF*) 

485  end 

H^EXP=[]; 

for  count=l:lengtii(plotw) 
if  meters 

waitbar(count/lengtli(plotw)); 

490  end 

2-^xpFk_exp+j*plotw(comit)*c_exp-plotvv(count)'^2*ni_exp; 

h_exp=mv(z__exp); 

h_exp=h_exp(aset,aset); 

skyindex=l; 

495  forindexl=l:aset_size 

for  index2=index  1  :aset_size 
redJioider(skyindex)=h_exp(mdexl,index2); 
skyindex=skyindex+ 1 ; 
end 

500  end 

H_EXP^P_,EXP  redjiolder']; 
end 

if  meters 

close(waitbarJiandle) 

505  end 

figure(5); 

plot(plotw/(2*pi),loglO(abs(H__EXP(ndx3d([l  aset_size*(aset_size+l )/2  length(plotw)],... 

510  1  ,skyred(cset_rel(mid_mdex),cset_rel(mid_index)),’:’)))),’r-’,. . . 

plotw/(2*pi),  Iogl0(abs(intplotl(ndx3d([aset_size  aset_size 
lengtli(plotw)],cset__rel(mid_index),cset_rel(mid_mdex),':’)))),’g-V.. 

plotw/(2*pi),  Iogl0(abs(intplot2(ndx3d([aset_si2e  aset_size 
lengtli(plotw)],cset_rel(mid__mdex),cset__rel(mid_index),*:')))),’b:*,. . . 

51 5  plotw/(2*pi),  Iogl0(abs(mtplot3(ndx3d([aset_size  aset__size 

length(plotw)],cset_rel(mid_index),cset_rel(mid_index),':')))),’k-.') 
ylabel([TI(*,int2str(cset(nnd_mdex)),V,mt2str(cset(inid_mdex)),*)  in/lbf (loglO  of)’]) 

xlabel('Omega  (Hz)’) 

520  grid  on 

"/otitieCEXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  1 INT  SOLUTIONS*) 
liii=legend(*Exp  FRF’, ... 

[num2str(lenctrl(l)),'  Hz  Bandwidtli  Corrected  FRF*],... 

[nmn2str(lenctrl(2)),'  Hz  Bandwidtli  Corrected  FRF],... 

525  [iium2str(lenctrl(3)X’  Hz  Bandwidtli  Corrected  FRF*]); 

axes(lili) 

print  -dmfile  fig5_5 
odlengtli=2*pi; 

530  countctrl=[3  10  50  200]; 

for  ind=l  :lengtli(countctrl) 
odivisions=countctrl(ind); 

owl  =odfreq-.  5*odlengtli:odlengtli/(odi  visions- 1  ):odffeq+.  5  ♦odlengtli; 

Wl=owl; 

535  %%compute  integrals  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%weiglit=ones(size(owl  ))*(omegax(  1 )); 
weiglit=owl; 
if  ind  ==  1 
owl 

540  df=(odlengtli/(odi  visions- 1  ))/(2  *pi) 
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end 

%%  compute  integrals  %%%%%%%%%%%%%%%%%%<>/o%o/m%%%%%%%%%%%%%%%%%%%o/oO/o 
DZ_R=[]; 

DZ_I=a; 

545  for  count=  1  :length(owl ) 

z__anal_red=kstat/(j  *owl  (00011!))+^  *owl  (count)*mstat; 
li_anal_red=inv(z_aiial__red); 

z_exp=k_exp/(j  *owl  (count))f  c_exp/(j  *owl  (count))+-j  *owl  (count)*m_exp; 
li_exp=inv(z_exp); 

550  li_exp=h_exp(aset,aset); 

liacc=h_aiial__red(cset_rel,cset_rel); 

hxcc=h_exp(cset_rel,cset_rel); 

dz=inv(inv(inv(hacc)*(hacc-lixcc)*inv(hacc))-hacc); 

dz_r=real(dz); 

555  dz_i=imag(dz); 

DZ_R=PZ_R  dz_r(:)]; 

,  DZJ=PZJ  dzj(:)]; 

end 

ifind==  1 
560  DZ_R; 

DZ^t 

for  i=l:lengtli(owl) 

dzterapl  =DZ_R(ndx3d([  1  cset_size*(cset_size+l  )/2  lengtli(owl )],  1 ,1  :cset__size*(cset_size+ 1  )/2,i)); 

565  dztemp2=DZ__I(ndx3d([l  cset_size*(cset_size+l)/2  lengtli(owl)],l,l:cset_size*(cset_si2ef  l)/2,i)); 

skyindex=l; 
for  index  1 = 1  :cset_size 
for  index2=indexl:cset_size 
570  dzl  (iiidexl  ,mdex2)=dztemp  1  (sl^iiidex); 

dzl(index2,iiidexl  )=dztempl  (skyiiidex); 
dz2(mdexl  ,index2)=dztemp2(skyindex); 
dz2(mdex2,mdexl  )=dztemp2(skymdex); 
skyindex=skyindex+l ; 

575  end 

end 
dzl 
dz2 
end 
i80 

end 

[INTK INTM  INTC]=intsub(DZ JJDZ_R,owl  ,W1 , 1  ./(j  ♦owl  ),cset_size); 
ifind==l 

[INTK  ESTTM  INTC]=mdsub(DZ_I,DZ_R,owl,Wl,l./(j*owl),cset_size); 

;85  INTM 
INTK 
INTC 
end 

90  ydM(m!^/o%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%o/o%o/^o/^o/^o/^o/^o/^o/^o/^o/^o/^o/^o/^^^^ 

%Qm)Stat 
kcorrected^kstat; 
mcorrected=mstat; 
ccorrected=cstat; 

95  kcorrected(cset_rel,cset_rel)=kcorrected(cset__rel,cset_rel)+INTK; 
mcorrected(cset__refcset_rel)=mcorrected(cset_rel,cset_rel>f ESTTM; 
ccorrected(cset_rel,cset_rel)=ccoiTected(cset_rel,cset_rel)+INTC; 
for  i=  1  :lengtli(plotw) 

tempza=kcoiTected+j*plotw(i)*ccorrected-plotw(i)'^2*mcorrected; 

DO  tempha=mv(tempza); 
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odcoiTha=[odcorTha  tempha(:)]; 
end 

if  ind==l 

mtplotl=odcoirha; 

605  save  INTMl  INTM 
save  INTKl  INTK 
saveINTCl  ESTTC 
elseif  md=2 
intplot2=odcorrha; 

610  save  INTM2 INTM 

save  INTK2 INTX 
save  INTC2  INTC 
elseif  ind==3 
intpiot3=odcon:ha; 

615  save  INTM3  INTM 

save  INTTC3 INTX 
save  INTC3  INTC 
else 

mtpIot4=odcorrha; 

620  save  INTM4  INTM 

save  INTK4  INTK 
save  INTC4  INTC 
end 
end 
625 

if  meters 

waitbarJiandle=waitbar(0,’Compiiting  Experimental  FRF'); 
else 

dispCGetting  experimental  FRF') 

630  end 

H_EXP=n; 

for  count=l  :lengtli(plotw) 
if  meters 

waitbar(comit/lengtli(plotw)); 

635  end 

z_exp=k__exp+j*plotw(count)*c_exp-plotw(count)^2*m_exp; 

h_exp=inv(z_exp); 

h_exp=h_exp(aset,aset); 

skyindex=l; 

640  for  indexl=l:aset_size 

for  index2=indexl:aset__size 
red_holder<skyindex)=li_exp(indexl  ,index2); 
skyindex=skyindex+ 1 ; 
end 

645  end 

H__EXP=[H_EXP  red_holder^]; 
end 

if  meters 

close(waitbarJiandle) 

650  end 

figure(6); 

plot(plot\v/(2*pi),logl  0(abs(H_EXP(ndx3d([  1  aset_size*(aset_size+l  )/2  length(plotw)],. . . 

655  I,skyred(cset_rel(mid_index),cset_rel(midjndex)),':')))),'r-V-* 
plotw/(2*pi),  Iogl0(abs(intplotl(ndx3d([aset_size  aset_size 
lengtli(plotw)],cset_rel(mid_index),cset_rel(mid__index),V)))),’g-',-.. 

plotw/(2*pi),  Iogl0(abs(intplot2(ndx3d([aset_size  aset_size 
lengtli(plotw)],cset_rel(midJndex),cset_rel(mid_index),':’)))),'b:V.. 
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660  plotw/(2*pi),  Iogl0(abs(mtplot3(ndx3d([aset_size  aset_size 

leiigtIi(plotvv)],cset_rel(niid_index),cset_rel(mid_index),’;')))),'k-.’,,.. 

plotw/(2*pi),  Iogl0(abs(intplot4(ndx3d([aset_size  aset_size 
iengtli(plotw)],cset_rel(mid_index),cset_rel(niid_mdex),':’)))),’iu:') 
ylabel(['H(’,int2str(cset(mid_index)V,',iut2str(cset(mid__mdex)),')  iii/lbf (loglO  of)’]) 

xlabel(’Oinega  (Hz)’) 
grid  on 

%title(’EXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  1 INT  SOLUTIONS’) 
lili=slegend(’Exp  FRF, ... 

37 0  [nuin2str(coimtctrl(  1 )),'  Pts  Corrected  FRF’], . . . 

[nuni2str(counlctrl(2)),'  Pts  Corrected  FRF’], ... 

[iiuin2str(comitctrl(3)),’  Pts  Corrected  FRF], ... 

[nuiii2str(countctrl(4)),’  Pts  Corrected  FRF]); 
axes(lili); 

575  print  -dinfile  fig5_6 
figure(7); 

plot(plotw/(2*pi),loglO(abs(H_EXP(ndx3d([l  aset_size*(aset_size+l)/2  lengtli(plotvv)],.., 

1  ,skyred(cset_rel(mid_mdex),cset_rel(imd_index)),':’)))),'r-’,. . . 

>80  plotw/(2*pi),  logl  0(abs(plotl  (ndx3d([aset_size  aset_size  lengtli(pIotw)],cset_rel(mid_index),cset_reI(mid__index),’:')))),’g- 

plotw/(2*pi),  Iogl0(abs(mtplotl(ndx3d([aset_size  aset_size 
Iengtli(plotw)],cset_rel(niid_mdex),cset__rel(inid_index),’;')))),’b:’) 
ylabel([’H(',int2str(cset(niid__index)),V,iiit2str(cset(niid_index)),’)  in/lbf (loglO  of)']) 

xlabel(’Omega  (Hz)’) 
if  titles 

%  title(EXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  1  OD&INTGRL  SOLTNS') 
end 

>90  grid  on 

lili=legend(’Exp  FRF, ... 

[nuiii2str(lenctrl(  1 )),'  Hz  Bandwidtli  MAT  Corrected  FRF’], ... 

{nuin2str(lenctrl(l)),'  Hz  Bandwidtli  INT  Corrected  FRF]); 
axes(lili) 

!95  if  pswitch=- y* 

print  -dcdjcolor 
end 

print  -dinfile  fig5_7 
00  figiire(8); 

plot(plotw/(2*pi),logl0(abs(H_EXP(ndx3d([l  aset_size’^(aset_size+l  )/2  length(pIotw)],... 
l,slcyred(cset_rel(mid__index),cset_rel(niid_index)),V)))),’r-',... 

plotw/(2*pi),  Iogl0(abs(plot2(ndx3d([aset_size  aset_size  lengtli(plotw)],cset_rel(mid_mdex),cset_rel(mid_index),':’)))),’g-- 

05  ’,... 

plotw/(2*pi),  iogl0(abs(intplot2(ndx3d([aset_size  aset_size 
lengtii(plotw)],cset_rel(mid_index),cset_rel(mid_index),':')))),'b:') 
ylabei(  fH(’,int2str(cset(niid_mdex)),V,int2str(cset(niid_index)),')  in/lbf  (log  1 0  of)’] ) 

1 0  xlabel('Omega  (Hz)’) 

if  titles 

%title(EXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  1  OD&INTGRL  SOLTNS’) 
end 
grid  on 

15  lili=legend(ExpFRF,... 

[nuin2str(lenctrl(2)),’  Hz  Bandwidth  MAT  Corrected  FRF], ... 

[nmn2str(lenctrl(2)),’ Hz  Bandwidtli  INT  Corrected  FRF]); 
axes(hli) 
if  pswitcli=='y' 
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720  print  -dcdjcolor 
end 

print  -dinfile  fig5_8 
figure(9); 

725  plot(plotw/(2*piXlog  1 0(abs(H_EXP(ndx3d( [  1  aset_size*(aset_size+ 1  )/2  lengtli(plot w)], . . . 

1  ,skyred(cset_rel(mid_index),cset_rel(niid_mdex)),V)))),V“V .  • 

plotw/(2*pi),  iogl0(abs(plot3(ndx3d([aset_size  aset_size  lengtii(plot\v)],cset_rel(mid_index),cset_rel(inid_index),V))))/g— 

5'  •' 

7 30  piotw/(2*pi),  log  1 0(abs(intpiot3(ndx3d(  [aset_size  aset_size 

lengtli(pIotw)],cset_relrinid_index),cset_rel(mid_mdex),':’)))),'b:’) 
ylabel(['H(',int2str(cset(inid_index)),’,',int2str(cset(mid_index));)  iii/lbf  (loglO  of)']) 

xlabel(’Omega  (Hz)’) 

735  if  titles 

rotiUeCEXPERMENTAL  FRF  AND  CORRECTED  FRFS  USING  MODE  1  OD&INTGRL  SOLTNS’) 
end 
grid  on 

Wi=legend(Exp  FRF', ... 

740  [nuin2str(lenctrl(3)),’  Hz  Bandwidtli  MAT  Corrected  FRF'], ... 

[nuin2str(lencfrl(3)),'  Hz  Bandwidtli  INT  Corrected  FRF']); 
axes(hli) 
if  pswitch==y 
print  -dcdjcolor 
745  end 

print  -dinfile  fig5_9 
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clear 

odivisioiis=3 
startloop=l 
skiploop=l 
enclloop=4 
cd  incomp  1 

save  odconf  odivisioiis  startloop  skiploop  endloop 

chap6_I 

cliap6_2 

cd .. 

clear 

odivisions=l 
startloop=3 
skiploop=l 
endloop=4 
cd  compl 

save  odconf  odivisioiis  startloop  skiploop  endloop 

chap6_l 

chap6__3 


119 


CHAP6.M 


20  %%%CHAP6_1%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%%% 

%%%Decompose  DZ  into  DM,  DK,  &  DC%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%usmg  matrix  formulation%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
25  %%%%%%%%%%%% 

clc 
clear 
closeall 
ioadINT 

30  %if  pswitch==y 
wlutebg(’wdTite’) 
close 
%end 
h=0; 

35  fignum=l; 

use__antires=l; 

titles=0; 

inid_index=round(lengtli(csety2); 

40  clear  L  Z_ANAL  c_exp  conn  h__anal_red 
clear  kexto  m_exp  mexto  temp_l_diags 

clear  true_dainping  true_mass  true_stiflhess  z_anal  z__anal_red  z_exp 
clear  L  Z_ANAL  c_exp  conn  h_anal_red 
clear  kexta  kexto  mex^  temp J_diags 
45 

%%%%%%%%%%%%%%%%%%%FORCE  THE  CSET%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%SPATIALLY1NC0MPLETE%%%%%%%%%%%%%%%%%% 
if  complete 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

50  %%%%%%%%%%%%%%%SPATIALLYCOMPLETE%%%%%%%%%%%%%%%%% 

cset=[4  5  6  7  8  9  10  11]%%%%%%%%%%%%%%%%%%%%% 
cseLrel=cset%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
cseLsize=lengtli(cset)%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%W/o%%%%%%%%%% 

55  else 

%%%%%%%%%%%%SPATrALLYINCOMPLETE%%%%%%%%%%%%%%%%%% 
cset=[l  3  5  7  9  11  13  15]%%%%%%%%%%%%%%%%%%%%%% 
cset_rel=[l  2  3  4  5  6  7  8]%%%%%%%%%%%%%%%%%%%%% 
cset__size=lengtli(cset)%%%%%%%%%%%%%%%%%%%%%%%%% 

60  end 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%®/o%>®/o%%%%%% 


load  exp 

65  odivisions=3 ; 

load  odcoiif 
nmnmodes=8  ; 
odlengtli=2.00 1  *pi ; 
firsttime=l; 

70  dw=odlengtii/(odivisions+l ); 

dovv=(freqtop-freqbottom)/(fmeness-l ); 
owref=freqbottom:dow:freqtop; 

if  firsttime 
75  if  meters 

waitbar  Jiandle=waitbar(0,'Computing  Experimental  FRF'); 
else 

dispCGetting  experimental  FRF') 
end 
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80  H_T^EXP-0; 

for  couiit=  1  :length(o\vref) 
if  meters 

waitbar(couiit/Iengtli(owref)); 

eiid 

85  z_exp=k_exp+j*o\vref(count)*c__exp-owref(couiit)'^2*m_exp; 
h_exp=mv(z_exp); 
h__exp=h_exp(aset,aset); 

%  skyiiidex=l; 

%  for  mdexl=l:aset_size 
90  %  formdex2=mdexl:aset_si2e 

%  red Jiolder(skyindex)=h_exp(iiidex  1  ,index2 ); 

%  skymdex=skyiiidex+l; 

%  end 
%  end 

95  %  H_T_EXP=pi_T^EXPredJiolder’]; 

H^T__EXP=p_T_EXP  li_exp(:)]; 
end 

if  meters 

close(waitbar_handle) 

100  end 

saveJH[_T_EXP  H_T^EXP 
clear  H_T^EXP 
%loadH_ANAL_RED; 

%iincorrha=H_ANAL_RED; 

105  %clearH^ANAL^RED 

uncorrha=[]; 
if  meters 

waitbarJiandle=waitbar(0,'Computing  Uncorrected  FRF'); 
else 

1 10  disp(’Getting  Uncorrected  FRF') 

end 

load  stat 

for  i=l  :lengtli(owTef) 
if  meters 

1 5  waitbar(i/lengtii(o’wref)); 

end 

tempza=kstat+j  *owef(i)*cstat-o\vref(i)^2*mstat; 

tempha=inv(tempza); 

uncorrha=[uncorrlia  teniplia(:)]; 

20  end 

if  meters 

close(waitbar__handle) 

end 

save  uiicorrha  uncorrha 
25  clear  uncorrha 

end 

loadH__T_EXP 

temp=[]; 

for  mdex=  1  :nummodes 
30  for  coord=l:aset_size 

if  index  ==  1 

antiresfreqs=fmd(owref>0  &  owref<omegax(index)); 
else 

antiresfreqs=fmd(owref>omegax(index-l)  &  ovvref<omegax(index)); 

35  end 

antires=mm(loglO(abs(H_T_EXP(udx3d([aset_size  aset^size  lengtli(ovvref)],... 
coord,coord,aiitiresfreqs))))); 

wlire=find(log  1 0(abs(H_T_EXP(ndx3d(  [aset^size  asel_size  lengtli(owref)] , . . . 
coord, coord,antiresffeqs))))==antires); 
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140  teiiip=ovvref(antiresireqs) ; 

aiitiresfreq(coord,index)=temp(wlire(  1 )); 
end 
end 

lemp=[]; 

1 45  for  mdex=  1  inummodes 

resfreqs=flnd(owref>omegax(index)-dow  &  owref<omegax(index)+dow); 
res=max(log  1 0(abs(H__T_EXP(ndx3d([aset_size  aset_size  length(o\vref)],. . . 
cset_rel(  1  ),cset_rel(  1  ),resfreqs))))); 

\vlire=find(logl0(abs(H_T_EXP(ndx3d([aset__sizeaset_size  lenglli(ovvreO],... 
1 50  cset_rel(  1  ),cset__rel(  I  ),resfreqs))))==res); 

tenip=owref(resfreqs) ; 
resfreq(index)=temp(wlire(  1 )); 
end 

clear  H_T_EXP 
155 

%clear  omegax 
%omegax=resfreq; 

for  loopindex=startloop:skiploop:endloop 
160  ODZ=[]; 

temp=[]; 
rdiag=[]; 
o\v=[]; 

lo\\'vver=ovvref(  1  )-.5*odIengtli; 

1 65  uppper=owref(  1  )+.5*odlengtli ; 

%  if  odivisions  ==  1 
%  o\v=[ow  owref(  1 )] 

%  else 

%  o\v=[ow  lowwer+dvvrdwiuppper-dw]; 

170  %  end 

%  rdiag=[rdiag  ones(l,odivisions*cset_size).*owref(l)]; 
for  index=l:loopmdex 

175  %  lowwer=omegax(index)-,5*odlenglli ; 

%  uppper=omegax(index)^-.5*odlength ; 

iovvwer=resfreq(index>.5  *odlength ; 
upppep=resjBreq(mdex)+.5*odlengtli ; 

180  o\v=  [ow  lowwer+dw:dvv:uppper-dw] ; 

rdiag=[rdiag  ones(  1  ,odivisions*cset_size).*resfreq(index)]; 
if  use_antires 

for  coord=  1  :aset_size 

1 85  lowwer=antiresfreq(coord,index)-.  5  *odlengtli; 

uppper=antiresfreq(coord,mdex)+.  5  *odlengtli; 

size(lowwer+dw:dw:uppper-dw) 

ow=[ow  lowwer+-dw:dw:uppper-dw] 

rdiag=[rdiag  ones(  1  ,odivisions*cset__size).  *resfreq(index)]; 

190  end 

end 
end 

R=diag(rdiag); 
load  stat 
1 95  if  meters 

\vaitbarJiandle=waitbar(0,’Computmg  DK/DM/DC’); 
else 

dispCGetting  DK/DM/DC’) 
end 


122 


CHAP6.M 


for  couiit=  1  :iengtli(ow) 
if  meters 

waitbar(comit/lengtli(ow)); 

end 

z__aial_red=kstat-o\v(couiit)^2*mstat; 

h_anal_red=mv(z_anal_red); 

z_exp=k_exp+j*ow(count)*c_exp-ow(coiuit)^2*m_exp; 

h_exp=iiiv(z_exp); 

h_exp=h_exp(aset,aset); 

liacc=h_anal__red(cset_rel,cset_rel); 

lixcc=li_exp(cset_rel,cset_rel); 

dz=mv(inv(mv(liacc)*(hacc-lixcc)*inv(hacc))-hacc); 

ODZ=[ODZ;dz]; 

temp=[temp;  ey^cset_size)  -ow(count)"'2*eye(cset_size)  +j*o\v(comit)*eye(cset_size)]; 
end 

if  meters 

close(waitbar_handle) 

end 

DKDME>Cl=inv(temp'*inv(R)*temp)*temp’*mv(R)*ODZ; 

%DKDMDCl=temp\ODZ ; 

ODKl=DKDMDCl(l:cset_si2e,:) ; 

ODMl=DKDMDCl(cset__size+l  :2*cset_size,:) ; 
OIX)l=DKDMDCl(2*cset_size+l:3*cset_size,:) ; 

clear  hacc  lixcc  li__exp  z_exp  z_anal__red  h_anal_red  dz 
clear  temp  ODZ 


%%%%%%%%%%%% 
load  stat; 
odcorrha=[]; 

kstat(cset__rel,cset_rel)=kstat(cset_rel,cset__rel)+ODKl; 
mstat(cset__rel,cset_rel)=mstat(cset_rel,cset_rel)+ODM  1 ;; 
cstat(cset__rel,cset_rel)=cstat(cset_rekcset_rel>4^DC  1 ; 

if  meters 

waitbar__handle=waitbar(0,'Coniputing  Corrected  FRF); 
else 

disp(’Getting  DK/DM/DC’) 
end 

for  i=l:length(owref) 
if  meters 

vvaitbar(i/lengtli(owref)); 

end 

tempza=kstaH-j*owref(i)*cstat-owref(i)^2*mstat; 
templia=inv(tempza); 
odcorrlia=[odcorrha  tempha(:)]; 
end 

if  meters 

close(waitbar_liandle) 

end 

if  loopindex  ==  1 
odcorrha  1  ^odcorrha; 
save  odl  odcorrhal 
clear  odcorrhal 
elseif  loopindex  ==  2 
odcorrha2=odcorTha; 
save  od2  odcorrha2 
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260  clear  odcoiTha2 

elseif  loopindex  ==  3 
odcorrha3=odcorTha; 
save  od3  odcorrha3 
clear  odcorrha3 

265  elseif  loopiiidex  ==  4 
odcorrha4=odcoiTlia; 
save  od4  odcorrha4 
clear  odcorrlia4 
elseif  loopiiidex  ==  5 
270  odcorrha5=odcorTha; 

save  od5  odcorrhaS 
clear  odcorrhaS 
end 

275  clear  tempha  tempza  kcorrected  mcorrected  ccorrected 

clear  ODMl  ODCl  DKDMDCl  ODZ  temp  kstat  mstat 

load  uncorrlia 
loadH^^T^XP 
280  if  complete 

plotwluchcoord=[cset(mid_index)]; 

else 

plotwliichcoord=[cset(mid_index)]; 

end 

285  for  plotindex=  1  :leiigtli(plotwliichcoord) 

which_coord=plotwliichcoord(plotmdex); 

li=figure(h+l); 

flgnum=fignum+l  ; 

subplot(2,l,l); 

290  plot(owrefi^(2*pi),.,. 

logl0(abs(H__T_EXP(ndx3d([aset__size  aset__si2e  lengtli(o\vref)]v 
wliich_coord,\viiich_coord,V)))),V-V . . 
owref/(2*pi),... 

Iogl0(abs(odcorrha(ndx3d([aset_size  aset_size  lengtli(owref)],... 

295  wliich_coord,  wliich_coord,':')))),'g“’,. . . 

ovvref7(2*pi),... 

Iogl0(abs(imcorrha(ndx3d([aset__size  aset_size  length(o\vref)],.- 
which_coord,wiiich_coord,’:')))),'b-.') 
ylabel(rH(’,int2str(aset(wliich_coord)),... 

300  V,uit2str(aset(which_coord)),’)  in/lbf  (loglO  of)']) 

hold  on 
v=axis; 

%xlabel(’C)mega  (Hz)') 
if  titles 

305  ifuse_antires 

title([int2str(odivisions);  PT  MATRIX  SOLTNS  WITH  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
else 

tiUe([mt2str(odivisions),'  PT  MATRIX  SOLTNS  WEIGHTED  LEFT  TO  RIGHT]) 
end 

310  end 

plot(resfreq(  1  :loopindex)/(2*pi),. . . 
ones(l,loopindex)*(v(3)+abs(v(4)-v(3))*.98),'k+') 

%  if  use_antires 

%  %plot(diag(antiresfreq(  1  rloopmdex,  1  :loopindex))/(2*pi),ones(  1  ,loopinde%x)*(v(3)+abs(v(4)-v(3))*.98),'k*’) 

315  %  end 

grid  on 

iiii=legend(’Experimental  FRF',’Corrected  MAT  Anal  FRF'/Uncorrected  Anal  FRF’,’hicluded  Modes'); 

axes(lili) 

hold  off 
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%if  pswitch==y 
%  prtfigl(fignuin) 

%  delete(h) 

%  %  else 

%  delete(h) 

%end 

subplot(2,l,2) 

dietop=fmd(abs(owref-ow(lengtli(ovv))*L5*ones(l,lengUi(owref)))==nmi(abs(owref- 

ow'ClengU^ow))*  1 . 5  *ones(  1  ,lengtli(owref))))); 
plot(owref(l  :tlietop)/(2*piX... 

loglO(abs(H_T_E)^(ndx3d([aset_size  aset_size  lengtli(owreO],... 
which_coord,wliich_coord,  1  :tiietop))));r-V.. 
owref(  1  :tlietop)/(2  *pi),. . . 

loglO(abs(odcorTha(ndx3d([aset__size  aset^size  lenglli(o\wef)],... 
wIuch_coord,whicli_coord,  1  :tlietop))))/g—V- 
owref(l  :thetop)/(2*pi),... 

Iogl0(abs(uiicorriia(ndx3d([aset_si2e  asetjsize  length(owref)],..- 
wliich_coord,wliich_coord,l  :tlietop)))),V/) 

ylabel([H(^mt2st^(aset(wluch_coord))/;4nt2st^(aset(wluch_coor^  iii/lbf  (loglO  of)']) 

hold  on 

v=axis; 

xlabel('Oinega  (Hz)’) 
if  titles 
if  use__aiitires 

title([int2str(odivisions);  PT  MATRIX  SOLTNS  WITH  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
else 

tiae([int2str(odivisions),'  PT  MATRIX  SOLTNS  WEIGHTED  LEFT  TO  RIGHT]) 
end 
end 

plot(resfreq(  1  :loopindex)/(2*pi),ones(  1  ,loopindex)*(v(3)+abs(v(4)-v(3  ))*.98),'k+’) 

%  if  use_antires 

%  %plot(antiresfreq(wliich_coord,l:loopindex)/(2*pi),ones(l,loopindex)*(v(%3)+abs(v(4)-v(3))*.98),’k*') 

%  end 

grid  on 

%lah=legend(’Experiniental  FRF’/Corrected  MAT  Anal  FRF’,TJncorrecled  Anal  FPr^Tncluded  Modes'); 

%axes(hli) 
hold  off 

if  loopindex  ==  2 
print  -dcdjcolor 
print  -dinfile  fig6_l 
elseif  loopindex  ==  3 
if  odivisions  ==  1 
print  -dcdjcolor 
if  complete 
print  -dmfile  fig6_9 
else 

print  -dmfile  fig6_6 
end 
end 

elseif  loopindex  ==  4 
print  -dcdjcolor 
print  -dmfile  fig6_2 
end 

%  if  pswitch=='y 

%  prtfigl(fignum) 

%  delete(h) 

%  else 
%  delete(h) 

%  end 
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380  end 

delete(h) 
clear  H^T^EXP 
clear  odcorrha 
clear  uncorrha 
385  end 
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%%%%GHAP6_2M%%%%%%%%0/oO/oO/o«/oO/oO/„0/,0/„0/„0/,0/„0/„0/„0/„0/„0/„0/„0/„0/„0/„0/„0/,^^^^^^ 

%%%%%%%%% 

%%%  Decompose  DZ  into  DM,  DK,  &  DC  %%%%%%V<P/o%%%%yoV«%%%%%%%%VoVo%%%V<PA%VoVoVo%%%% 
%%%  using  matrix  fonnulation  %%%yo%%%VoVo%VoVoVoVo%V<PmVoV<,%Vo%%%y,^^^^^^^ 

clc 

clear 

closeail 

loadINT 

whitebgCwhite') 

close 

h=0; 

figiiuiii=l; 

use_aiUires=l 

titles=0 

mid_mdex=roimd(lengtli(cset)/2); 

clear  L  Z_ANAL  c_exp  couii  h_anal_red 
clear  kexto  m_exp  mexto  temp J^diags 

clear  true^damping  true_mass  true^stiffness  z__aiial  z_aiial_red  z_exp 
clear  L  Z_ANAL  c__exp  comi  h_anal_red 
clear  kexta  kexto  mextO  temp_l_diags 

%VoVoVo%Vo%%Vo%%VoVo%%%VoV<p/omRC'ETmcsE^^ 

%%%%%%%%%%%%SPATIALLYINC01^LETE%%%%%%%%%%o/oO/oO/oO/^o^o/^o^o/^ 

if  complete 

VoVoVoV<P/o%Vo%VoViP/o%Vo%%%VoViP/oViP/o%%%VoVoVo%VoVo^ 

%%%%%%%%%%%%%%%SPAT[ALLY  COMPLETE%%%%%%%%%%%%%%%o/oO/o 
cset=[4  5  6  7  8  9  10  ll]%%%%%%%%%%%%%%%%%%%o^% 
cset_^el=cset%%%%%%%%%%%%%%%%%%%%o/o%%%o^%o^%%%%%o/^^o^^ 
cset_size=length(cset)%%%%%%%%%%%%%%%%%%%%%o/oO/^o^o/^o/^o^ 
%%%%%%%%%%%%%%%%%%%%%%0/o0/o%%%0/o%0/o0/o%%%0/o%%0/o%%0/o%%%0/o%0/o%0/o0/o 

else 

%%%%%%%%%%%o/oSPATL^LYINCOMPLETE%%%%%%%%%%%%%%%%%o/o 
cset=[l  3  5  7  9  11  13  15]%%%%%%%%%%%%%%%%%%%%%% 
cset_rel=[l  2  3  4  5  6  7  8]%%%%%%%%%%%%%%%%o/o%%%% 
cset_size=lengtli(cset)%%%%%%%%%%%%%%%%%%%%%o/oO/o%o/o 
end 

%%%%%%%%%%%%%%%%%%%%%%%%%%o^o^o^%%o^o^%o^o^o^%y^^%o^%%o^%o^o/^,o^o^y^ 

%%%%%%%%%%%%%%%%%%%%0/oO/oO/oO/oO/^0/^0^0^0/^0^0/^0/^0/^0/^0/^^^^^ 

%%%%%%%%%%%%%%%SPATIALLYCOMPLETE%%%%%%%%%%%%%%%%0/0 
%cset=:[7  8  9  10  11  12  13  14]%%%%%%%%%%%%%%%%%o/oO/oO/oO/o 
%cset_^rel=cset%%%%%%%%%%%%%%%%%%%%%o/o%%%%%o/o%%%%%%o/o% 
%cset^size=lengtli(cset)%%%%%%%%%%%%%%%%o/o%%%%%%o/oO/o%% 

%%%%%%%%%%%%%%%%%%%%%%0/o%%0/o%%%0/o%0/o%%%0/o%0/o%%0/o0/o%%0/o%0/o%%0/o 

load  exp 
odivisions=3  ; 
load  odconf 
niiiimiodes=9  ; 
odlengtIi=2*pi ; 
firsttime=l; 

dw=odlengtli/(odivisions+l ); 
dow=(freqtop-freqbottom)/(fmeness-l ); 
owref=freqbottom:dow:freqtop; 

if  firsttime 
if  meters 
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waitbarJiandle=waitbar(0;Computing  Experimental  FRF*); 
else 

disp('Getting  experimental  FRF') 

450  end 

H^T_EXP=[]; 
for  comit=l  :lengtli(owref) 
if  meters 

455  waitbar(count/lengtli(ovvref)); 

end 

z_exp=k_exp+j*owref(coimt)*c_exp”Ovvref(count)^2*m_exp; 

li_exp=inv(z_exp); 

h_exp=h_exp(aset,aset); 

460  %  skyindex=l; 

%  for  index  1 = 1  :aset_size 

%  for  index2=indexl  :aset_size 

%  red Jiolder(skyindex)=h__exp(index  1  ,index2); 

%  skyindex=skyindex+l; 

465  %  end 

%  end 

%  H_T^EXP=[H_,T_EXP  red Jiolder']; 

H^T^EXP=[H__T_EXP  h^exp(:)]; 
end 
470 

if  meters 

close(waitbarJiandle) 

end 

475 

save  H^T^EXP  H_T_EXP 
clear  H^T^XP 
end 

480  %loadH^ANAL^RED; 

%uncorrha=H_ANAL_RED; 

%clear  H_ANAL_RED 
if  firsttime 
micorrha=[]; 

485  if  meters 

waitbar__handle=waitbar(0;Computing  Uncorrected  FRF'); 
else 

dispCGetting  Uncorrected  FRF') 
end 
490 

load  stat 

for  i=  1  :length(owref) 
if  meters 

waitbar(i/lengtli(owref)); 

495  end 

tempza=kstat+j  *owref(i)*cstat-owref(i)^2  *mstat; 
tempha=mv(tempza); 
uncorrha=[micorrha  templia(:)]; 
end 
500 

if  meters 

close(waitbar_handle) 

end 

505  save  uncorrha  uncoirha 
clear  imcorrha 
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end 

loadH_T_EXP 

510  for  mdex=  1  rnuminodes 

for  coord=  1  :aset_size 
if  index  ==  1 

antiresffeqs=find(owref>0  &  owref<omegax(index)); 
else 

5 1 5  antiresffeqs=rmd(owref>omegax(index- 1 )  &  owref<omegax(index)); 

end 

aiitires=niiii(logl0(abs(H_^T_.EXP(ndx3d([l  aseLsize*(aset_size+l )/2  length(ovvreO],... 

1  ,skyred(coord,coordXantiresfreqs))))); 

wIire=fmd(logl0(abs(H_T_,EXP(ndx3d([l  aset_size*(aset_sizefl )/2  lengUT(o\vref)],... 

520  1  ,skyred(coord,coord),antiresfreqs))))==antires); 

temp=owTef(antiresffeqs) ; 
antiresfreq(coord,index)=temp(whre); 
end 
end 
525 

for  indeX“l:iiuimnodes 

resfreqs=fmd(owref>oniegax(index)-dow  &  owref<omegax(index>f dow); 
res=max(log  1 0(abs(H__T_EXP(ndx3d([aset__size  aset__size  length(o\vref)],. , . 
cset_rel(  1  ),cset j-el(  1  )^esif eqs))))); 

530  wlire=find(logl0(abs(H_T_EXP(ndx3d([aset_size  aset__size  length(owreO],-. 
cset_rel(  1  ),cset_rel(  i  ),resfr eqs))))==res); 
teinp=ovvref([resfreqs) ; 
resffeq(index)=temp(wlire); 
end 
535 

clear  H_T_EXP 
%clear  omegax 
%omegax=resffeq; 

540  for  loopindex=startloop:skiploop:endloop 

subdivisions=9;  %  should  be  odd  for  sinipson’s  rule 
intlengtli=2*pi;  %  1  Hz  bandwidtii 
dW=intlengtli/subdivisions;  %  sampling  freq  =.25  HZ 
W=[]; 

545  Wl=[]; 

vveiglit=[]; 

lowvvei=o\vref(  1  )-,5*odlengtli; 
uppper=owref(  1  >f.5*odlength; 

W=[W  lo\vwer4-dW:dW:uppper-dW]; 

550  W1=[W1  lowvverfdW:dW:uppper-dW]; 

weiglit=[weiglit  ones(size(lowweiH-dW:dW:uppper-dW))]; 

for  index=  1  doopindex 
low\ver=omegax(index)-.5*intlengtli; 

)55  upppei=omegax(index)+.5*intlengtli; 

W=[W  lowwerfdW:dW:uppper-<iW]; 

W1=[W1  lowwei+dW:dW:uppper-dW]; 
weight=[weiglitones(size(lowweri-dW:dW:uppper-dW))]; 
lo\ver=antiresffeq(index)--5*intlengtli; 
i60  upper=antiresfreq(index)+.5*intlength; 

W=[W  lowwerfdW:dW:uppper-dW]; 

W1=[W1  lowwerJ-dWidWiuppper-dW]; 
weiglit=[weightones(size(lowwer+dW:dW:uppper-dW))]; 
end 
65 

%%cx>mputemtegrals%%%%%%%%%%%%%%%%%%%o/oO/oO/o0^o^o^o/^o/^o^o/^o/^o/^o/^o/^o/^o/^o/Qo^o^^ 
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DZ_R-[]; 

DZJ=[]; 
if  meters 

570  waitbar Jiaiidle=waitbar(0;Computiiig  DZ'); 

else 

disp(’Computing  DZ'); 
end 

load  stat 

575  for  coiint=  1  :lengtli(  W) 

if  meters 

vvaitbar(count/length(W)); 

end 

z_^anal_red=kstat/(j*W(comit))+cstat/(j*W(count)>+-j*W(coimt)*mstat; 
580  li_anal_red=inv(z_anal__red); 

hacc=h_anal__red(cset_rel,cset_rel); 

z_exp=k_exp/ (j  *  W(coiuit))+-c_exp/(j  *  W(comit))+j  *  W(count)*m_exp; 
h__exp=inv(z_exp); 

585  h__exp=h_exp(aset,aset); 

lixcc=h__exp(cset_rel,cset_rel); 

dz=inv(inv(iiiv(hacc)*(hacc-lL\cc)*inv(hacc))-hacc); 

dz_r=real(dz); 

590  dzj=imag(dz); 

DZ^R=PZ^R  dz_r(:)]; 

DZJ=pZJdzJ(:)]; 

end 

595  if  meters 

close(waitbarJiandle) 

end 

W1=1.AV1; 

[DSTTK,  INTM,  INTC]=iiitsub(DZ^I,DZ_R,W,Wl,weiglit,cset__size); 

600  intcorrha=[]; 

clear  hacc  lixcc  li_exp  z__exp  z_anal_red  h_anal__red  dz 
clear  temp  ODZ 


605  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%o/o%%%%o/o0/o%%%%o/o%%%%%% 

%%%%%%%%%%%% 
load  stat; 
odcorrlia=[]; 

610  kstat(cset__rel,cset_rel)=kstat(cset_rel,cset__rel)+INTK; 

mstat(cset_rel,cset_rel)=mstal(cset_rel,cset_rel>i’INTM; 

cstat(cset_rel,cset_rel)=cstat(cset__rel,cset__rel>i-INTC; 

if  meters 

615  \vaitbarJiandle=waitbar(0, 'Computing  Corrected  FRF'); 

else 

disp(’Getting  DK/DM/DC') 
end 

620  for  i=  1  ;lengtli(owref) 

if  meters 

waitbar(i/lengtli(o\vref)); 

end 

tempza=kstat+j  *owref(  i  )*  cstat-owref(i  )'^2*mstat; 

625  templia=inv(tempza); 

odcorrha=[odcoiTha  tempha(;)]; 
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end 

if  meters 

close(waitbarJiaiidle) 

end 


if  loopindex==  1 
lodcorrha  1  =odcoiTha; 
save  lodl  Icxlcorrhal 
clear  lodcorrha  1 
elseif  loopindex  ==  2 
lodcorriia2=odcorrha; 
save  lod2  lodcorrha! 
clear  lodcorrha! 
elseif  loopindex  ==  3 
lodcorrha3=odcoiTha; 
save  Iod3  lodcorrhaS 
clear  lodcorrliaS 
elseif  loopindex  ==  4 
lodcorrha4=odconha; 
save  lod4  lodcorrha4 
clear  lodcorrha4 
elseif  loop'index  ==  5 
lodcorrha5=odcorrha; 
save  lod5  lodcorrhaS 
clear  lodcorrliaS 
end 


clear  tempha  tempza  kcorrected  nicorrected  ccorrected 
clear  MIM INTC INTIC  temp  kstat  mstat  cstat 
load  micoirha 
loadH__T_EXP 

if  complete 

plotwliichcoord=[cset(mid_mdex)]; 

else 

plot\\iiichcoord=[cset(mid_index)] ; 
end 

for  plotindex=l:length(plotwliichcoord) 
wliich_coord=plotwhichcoord(plotindex); 
h=figure(h+l); 
subplot(!,lJ) 
plot(owref7(2*pi),.. , 

loglO(abs(H_T_EXP(ndx3d([aset_size  aset_size  lengtli(owref)],.-. 
\vliich_coord,wliich_coord,V))));r-',. . . 
owref/(!*piX... 

Iogl0(abs(odcorrha(ndx3d([aset__size  aset^size  lengtli(owref)],.-. 
wIiich_coord,which_coord,':’))));g-V  .. 

owref7(!*pi),... 

loglO(abs(micorTha(ndx3d([aset_si2e  aset^size  length(owreO],... 
wliich_coord,wlnch_coord;:')))),'b-.') 
ylabel([H(',int2str(aset(>vluch_coord)),. . . 

V,iut!str(aset(wliich_coord));)  in/lbf (loglO  of)’]) 
hold  on 
v=axis; 

%xlabel(’Omega  (Hz)') 
if  titles 

if  use_antires  ==  'on ' 
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title([mt2stii(odivisions),'  PT  INTEGRAL  SOLTNS  WTl'H  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
else 

title([mt2str(odivisions);  PT  INTEGRAL  SOLTNS  WEIGHTED  LEFT  TO  RIGHT!) 

690  end 

end 

plot(omegax(  1  :loopindex)/(2  *pi),. . . 
ones(  1  ,loopmdex)*(v(3)4-abs(v(4>v(3))^.98),'k+') 

%  if  use__antires 

695  %  %plot(antiresfreq(wliich_coord,l:loopindex)/(2*pi),ones(l,loopindex)*(v(%3>fabs(v(4>v(3))*.98),'k*’) 

%  end 
grid  on 

hh=iegend('Experimental  FRFVCorrected  INT  Anal  FRF',TJncoiTected  Anal  FRF^Tncliided  Modes’); 
axes(hli) 

700  hold  off 

%if  pswitch==’y‘ 

%  prtfigl(fignuni) 

%  deiete(h) 

%else 

705  %  delete(h) 

%end 

subplot(2,l,2) 

thetop=fmd(abs(owref-W(length(W))*  1 .3*ones(l,length(owref)))==niin(abs(owref- 
W(lengtli(W))*  1 .3*ones(l,lengtli(owre0)))); 

710  %  h=figure(h+l ); 

%  fignmn=fignmn+l  ; 
plot(owref(  1  :tlietop)/(2*pi),. . . 

loglO(abs(H_T_E^(ndx3d([aset_size  aset_size  length(owref)],.-- 
wliich_coord,wliich_coord,  1  :tiietop)))),'r-’,. . . 

715  owref(l:tlietop)/(2*pi),... 

loglO(abs(odcoiTha(ndx3d([aset_size  aset_size  lengtli(o\\Tef)],... 
wliich^coord,  wliich_coord,  1 :  tlietop)))),’g~’, . . . 
owrefi;i  :thetop)/(2*pi),... 

logl0(abs(xmcoiTiia(ndx3d([aset_size  aset_size  lengtli(o\vref)],... 

720  which_coord,wliich_coord,  1  :tlietop)))),’b-,’) 

ylabel([TI(',int2str(aset(wliich_coord)),’,',int2str(aset(which_coord)),’)  in/lbf  (log  1 0  oi)'] ) 

hold  on 

v=axis; 

xlabel(’Omega  (Hz)’) 

725  if  titles 

if  use_antires  ==  'on ' 

title([int2str(odivisions);  PT  INTEGRAL  SOLTNS  WITH  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
else 

tiae([int2str(odivisions);  PT  INTEGRAL  SOLTNS  WEIGHTED  LEFT  TO  RIGHT]) 

730  end 

end 

plot(omegax(l:loopindex)/(2*pi),ones(Uloopindex)*(v(3>t-abs(v(4>v(3))*.98X’k+') 

%  if  use__antires 

%  %plot(antiresfreq(which_coord,  1  :loopindex)/(2*pi),ones(  1  ,loopindex)*(v(%3)f  abs(v(4)-v(3))*.98),'k*’) 

735  %  end 

grid  on 

%lili=legend(’Experiinental  FRF',’Corrected  INT  Anal  %RF','Uncorrected  Anal  FRF’,'hicluded  Modes’); 

%axes(lili) 

740  if  loopindex  ==  2 

print  -dcdjcolor 
print  -dinfile  fig6_3 
h=figure(h+l); 
hold  off 
745  load  od2 

plot(ovvref(  1  ;tlietop)/(2  *pi), . . . 
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loglO(abs(H_T_EXP(ndx3d([aset_si2e  aset_size  length(owref)],... 
wluch_coord,which_coord,l:tlietop)))),'r-',... 

owref(l  :tlietopy(2*pi),... 

Iogl0(abs(odcorrha(ndx3d([aset_size  aset_size  lengtli(o\vre£)],... 

wiiich_coord,\vliich_coord,l;tlietop)))),'g-',... 

owref(l  :aietopy(2*pi),... 

Iogl0(abs(odcorrha2(ndx3d([aset_size  aset_size  leiigtiyowreO],- 
which_coord,which_coord,l  itlietop)))),^.’,... 
owref(l  ;thetopy(2*pi),.  .. 

Iogl0(abs(uncorrha(ndx3d([aset_size  aset_size  lengtli(owreO],- 
which_coord,which_coord,  1  :aietop)))),'m:') 

ylabel([H(^mQstr(aset(w1iich_coord));;,mt2str(aset(\vliich_coord)),')  in/lbf (loglO  of)']) 

hold  on 
v=axis; 

xlabel('Oinega  (Hz)') 
if  titles 

if  use_antires  =  'on ' 

tiUe([int2str(odivisions),'  PT  INTEGRAL  SOLTOS  WITH  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 

title([int2str(odivisions),'  PT  INTEGRAL  SOLTNS  WEIGHTED  LEFT  TO  RIGHT]) 
end 
end 

pIot(omegax(  1  :loopindex)/(2  *piXones(  1  ,loopindex)*(v(3>f abs(v(4>v(3))*,  98),'k+') 

%  if  use_antires 

%  ®^Plot(aittiresfreq(wMch_coord,l:loopmdexy(2*pi),ones(l,loopindex)*(v(%3)fabs(v(4>v(3))*.98),'k*') 

%  end 

grid  on 

hli=legend(’Experimental  FRFVCorrected  INT  Anal  FRF/Coirected  MAT  Anal  FRF'.TJncoiTected  Anal  FRF  Included 
Modes’); 

axes(lih) 
hold  off 
pmit  -dcdjcolor 
print  -dinfile  fig6_5 
elseif  loopindex  =  3 
if  odivisions  ==  1 
print  -dcdjcolor 
print  -dinfile  fig6_7 
h=figure(h+l); 
hold  off 
loadod3 

piot(owref(l  :tlietop)/(2*pi),... 

loglO(abs(H_T_]E^(ndx3d([aset_size  aset__size  length(owref)],... 

\vliich_coord,which_coord,l:tlietop)))),'r-’,.,. 
owTef(  1  :tlietop)/(2*pi),... 

Iogl0(abs(odcorrha(iidx3d([aset_size  aset_size  Iengtli(ov\Tef)],... 

\riiich_coord,wluch_coord,  1  ;tlietop)))),’g-*V 
owref(  1  :tiietop)/(2*pi),.,. 

Iogl0(abs(odcorrha3(ndx3d([aset_size  aset_size  lengtli(o\vref)]v  • 
which_coord,\vlnch_coord,  1  itlietop)))),!)-.’,. . . 
owref(  1  :tlietop)/(2*pi),. . . 

loglO(abs(uncorTha(ndx3d([aset_size  aset^size  lengtli(owref)],... 

which_coord,wliicli_coord,  1  itlietop)))),^:’) 

ylabel(['HC,iiit2str(aset(wluch_coord));;,mt2str(aset(\vhich_coord));)  in/lbf  (log  10  ot)’]) 

hold  on 

v=axis; 

xlabel(’Omega  (Hz)’) 
if  titles 

ifuse__antires  ==  ’on’ 

^tiae([mt2str(odivisions),'  PT  INTEGRAL  SOLTNS  WITH  ANTIRESONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
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820 


825 


830 


835 


title([int2str(odivisions);  PT  INTEGRAL  SOLTNS  WEIGIilED  LEFT  TO  RIGHT]) 
end 
end 

plot(omegax(  I  :loopmdex)/(2*piXones(  1  ,loopindex)*(v(3)+abs(v(4)-v(3))*.98);k+') 

%  if  iise_antires 

%  ^®Plot(antiresfreq(wliich_coord,l:loopindex)/(2*pi),ones(l,loopindex)*(v(%3>4-abs(v(4)-v(3))*,98),’k*') 

%  end 

grid  on 

iili^legendCExperiniental  FRF’;Corrected  INT  Anal  FRFVUncorrected  Aiial  FRl^TncIiided  Modes')- 

axes(lih) 

hold  off 

print  -dcdjcolor 

print  -dmfile  fig6_8 

end  ^ 

elseif  loopindex  ==  4 
print  -dcdjcolor 
print  -dmfile  fig6_4 
end 

%  ifpswitch=='y' 

%  prtfigl(fignum) 

%  delete(h) 

%  else 
%  delete(h) 

%  end 
end 

delete(h) 
clear  H^T_EXP 
clear  odcorrha 
clear  uncorrha 
end 


loadH_T_EXP 
840  load  odl 

loadod2 
loadod3 
load  od4 
load  uncorrha 
845  %load  lod5 

h=figure(h+l ); 
plot(owref7(2*pi),... 

loglO(abs(H_T_EXP(ndx3d([aset_size  aset_size  lengtli(owref)],... 
which_coord,wliich_coord,':’)))),'r- . . 

850  ovvref7(2*pi),.., 

Iogl0(abs(odcorrhal(ndx3d([aset_si2e  aset_size  lengtli(owref)],  • 
wliich_coord,which_coord,':')))),’g— , . 
owref/(2*pi),... 

Iogl0(abs(odcorrha2(iidx3d([aset_size  aset_size  lengtli(ov\Tef)],,.. 
855  which^coord,  which_coord,':')))),'b-.  V . . 

o\vref7(2*pi),... 

logl0(abs(odcoiTha3(ndx3d([aset_size  aset_size  length(o%^Tef)],... 
wliich_coord,wliich_coord,':')))),'m-', . . . 
o\vref7(2*pi),... 

860  Iogl0(abs(odcorrha4(iidx3d([aset_si2e  aset__size  length(o\vref)],-- 

wliich_coord,wiiich_coord,':')))),'c-.  V . . 
o\vref7(2*pi),... 

Iogl0(abs(uncorrha(ndx3d([aset_size  aset_size  lengtli(o\vref)],-.- 
which_coord,wliich_coord,':’)))),'k-.') 

865  ylabel([’H(',int2str(aset(whicli_coord)),... 

’,',int2str(aset(\vliich_coord)),’)  in/lbf  ( log  1 0  of)'] ) 
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hold  on 
v=axis; 

%xlabel(’Oniega  (Hz)') 

870  if  titles 

if  use_antires  ==  'on ' 

title([mastr(odivisions),’  PT  INTEGRAL  SOLTNS  WITH  ANTD^SONANCES  WEIGHTED  LEFT  TO  RIGHT]) 
else 

title([iiit2str(odivisioas);  PT  INTEGRAL  SOLTNS  WEIGHTED  LEFT  TO  RIGITTl) 

875  end 

end 

%  plot(omegax(  1  :loopindex)/(2  *pi ),. . . 

%  ones(  l,loopindex)*(v(3)+abs(v(4)-v(3))*.98),'k+') 

%  if  use_antires 

880  %  %plot(antiresfreq(which_coord,l:loopindex)/(2*pi),one$(l,loopindex)*(v(%3)+abs(v(4>v(3))*.98),'k*') 

%  end 
v=axis; 

axis([10  300  v(3)  v(4)]); 

'  grid  on 

385  lili=legendCExperimental  FRF','!  mode  MAT  solutiou’,'2  mode  MAT  solution’,3  mode  MAT  soIution';4  mode  MAT 

solution', TJncoixected  Anal  FRF'); 
axes(lili) 
hold  off 

%ifpswitch==’y’ 

190  %  prtfigl(fignum) 

%  delete(h) 

%else 

%  delete(h) 

%end 

195 

print  -dcdj  color 
print  -dmfile  fig6_12 

00 
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%  PROGRAM:  SETUP.M  o/, 

%  INimiZES  DATA  FILE  FOR  A  FINITE  ELEMENT  ANALYSIS 
%  VARIABLES  SAVED  TO  FILE  SETUP.MAT 


clc 

clear 

L=input(’Enter  tlie  total  lengdi  of  tlie  beani:  ’); 
nuiiiel=mput('Nuin  of  elements: '); 
dof_node=mput(’num  dof  per  node:  ’); 
area=mput(’area  of  beam:  ’); 
eeii=mput('EI  of  beam: '); 

%ee=input(’modulus  for  beam:  ’); 
plio=input(’mass  density  for  beam:  ’); 
nmndof^dof_node*(numel+ 1 ) 
clc 

wliile  1 

conforce=input('Enter  #  of  concentrated  loads: '); 
if  ~isnan(conforce) 
break 
end 
end  • 

for  i=l:conforce 

forcepos(i)=input([’Beam  position  of  force  ',nura2str(i),' : ']); 
forcesiz(i)=input(['Magnitude  offeree  ',num2str(i); :  ’]); 
end 

force=2eros(nmndof,  1 ); 
for  i=l:conforce 

temppos=ceil(forcepos(i)/(L/nmnel)); 

force(dof_node*temppo^l  )=force(dof_node*temppos+l  .5*forcesiz(i); 
force(dof_node*temppos+3)=force(dof_node*temppos+3>f.5*forcesiz(i); 

end 
while  1 

hnass=input('Enter  #  of  Limiped  masses: '); 
if--'isnan(lmass) 
break 
end 
end 

for  i=l:imass 

masspos(i)=mput(['Beam  position  of  Mass  ’,num2str(i),' : ']); 
masssi2(i)=input([’Magnitude  of  Mass  ’4ium2str(i),’ : ']); 
end 

lumpmass=zeros(nmnel,  1 ); 
for  i=l:lmass 

temppos=round(masspos(i)/(L/numel)); 

lmnpmass(temppos>=lumpmass(temppos)+masssiz(i); 

end 

lumpmass’ 
wliile  1 

lspring=input(’Enter  #  of  Lumped  Springs:  *); 
if~isnan(lspring) 
break 
end 
end 

for  i=  Lispring 

springpos(i)=input([’Beain  position  of  Spring  ’,num2str(i),’ :  ’]); 
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spriiigsiz(i)=mput(['Sprmg  constant  of  Spring  ',nuni2str(i); : ']); 
end 

iimipsprmg=zeros(nuniel,l ); 
fori=l:lsprmg 

65  temppos==rouiid(sprmgpos(i)/(L/nuinel)); 

lunipspring(temppos)=luinpspriiig(temppos>fspringsi2(i); 

end 

lumpspring’ 

cx>iin=[U]; 

70  fori=2;numel 

conn=[conn;i,i+l]; 
end 
comi 
while  1 

75  bc=[’pimied-pinned  ' 

'clamp-clamp  ' 
left  guided  clamp ' 

'riglit  guided  clamp’ 

'cantilevered  ' 

80  'free-free  ']; 

clc 

help  bctext 

n=input('Select  a  boundary  condition: '); 
if((n>0)&(n<7)) 

85  break 
end 
end 

ifn==  1 
bc='pp'; 

90  elseif  n  ==  2 

bc='cc’; 
elseif  n  =  3 
bc='ic'; 
elseif  n  ==4 
95  bc='rc’; 

elseif  n==5 
bc=’cl' 
else 
bc=’fF; 

)0  end 

clear  i 

clear  temppos 
clear  conferee 
)5  clear  forcepos 

clear  forcesiz 
clear  Imass 
clear  masspos 
clear  masssiz 
0  clear  Ispring 

clear  springpos 
clear  springsiz 
save  setup.mat 
x=0:Lynumel:L; 

5  for  i=l:length(x) 

li(i)=5; 
end 
clg 

hold  off 

3  %axis(’ofF) 
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%piot(n,[]) 
%hold  on 
plot(x4i,x,li,'x') 
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%%yo%%%%%%yoyo%%%%%%%%%%%%yoyo%%%%%yo%%%%%%%%%yo%%%%%%%%%%%%yo%%%%%% 

%%%%%%% 

%FORM  FE  BEAM  MODEL  WITH  ERRORS  AT  PRESCRIBED  LOCATIONS%%%%%%%%%%% 

%MASS  ERROR  LOCATION  SPECIFIED  VIA  POSOFMASSER.  VALUE  OF  ERROR  IN%% 

5  %PERCENT  GIVEN  BY  VALOFMASSERR,  ETC.%%%%%%%%%%%%%%%%%%o^%%%%%%%o^o/oO^%% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%% 

clear,clc;clg 

casenanie='Spatially  Incomplete  A  set  Stiffness  error'; 

10  ifexistCsetup.mat') 

load  setup  %load  existing  input  data  file 
else 

setup  %create  new  input  data  file 

end 

1 5  struc_damping=0.0000000000 1 ; 

Iunipdamp=luinpspring*0  ; 
posofMasseiT=[3  4]  ; 
valueofMasseiT=[0.25  0.25]  ; 
posofStifrerr=[3  4]; 

20  valueofStifrerr=[0.25  0.25]  ; 
posofDamperr=[3  4]; 
valueofDainperr=[0.25  0.25]  ; 
for  i=l:lengUi(posofMasserr) 
lunipniass(posofMasserr(i))=valueofMasseiT(i); 

25  end 

for  i=l:lengtli(posofStifrerr) 
luiupspring(posofStifferr(i))=valueofStifferr(i); 
end 

for  i=l:lengtli(posofDamperr) 

30  lmnpdanip(posofDamperr(i))=valueofDampeir(i); 

end 

elen=L/numel; 
ifbc=-pp' 
do£tokill=2;  ' 

35  elseifbc— cc' 

doftokill=4; 
elseifbc==’fr 
doftokill=0; 
elseifbc==’cr 
40  doftokill=2; 

end 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%0/o%0/o%%%%0/o%0/o%%0/o%%%%%0/o0/o 

%BUILD  THE  ELEMENTAL  MASS  AND  STDFFNESS  MATRICES%%%%%%% 

^5  %FOR  A  2DOF/NODE  "F.E."  STRUCTURES.%%%%%%%%%%%%%%%%%%%%% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%0/o0/o0/o0/^0/^0/^0/^0/^0^0/^0^0/^0/^0^^^^ 

%ELEMENT  STIFFNESS  MATm%%%%%%%%%%%%%%%%%%%%%%%%o/oO/oO/^o/^^o^o^o/^ 
ke=[12  6*elen-12  6*elen; 
iO  6*elen  4*(elen^2)  -6*elen  2*(elen^2); 

-12  -6*elen  12  -6*elen; 

6*elen  2*(elen'^2)  -6*elen  4*(elen^2)]; 

ke=(eeii/(elen'^3 )).  *ke; 
i5  g=386; 

elemke=eeii/(elen^3); 

^ELEMENT  MASS  MATmCE%%%%%%%%%%%%%%%%%%%%%%%%%0/0 
me=[156  22*elen  54  -13*elen; 

0  22*elen  4*(eien'^2)  13*elen  -3*elen^2; 
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me  =  (((pho*area)*elen/g)/(420)).*me; 

65  elemme=(((pho*area)*eleii/g)/(420)); 

%%%%%%%%%%%% 

%ASSEMBLE  GOBAL  STIFFNESS  AND  MASS  MATRIX  FOR  A  2  DOF  F.E.  STRUCTURE%%%% 

70  %BASED  ON  THE  ELEMENTAL  MATRIXES  KE  AND  ME.  THE  STRUCTURE  CONSISTS  OF%%% 

%NUMEL  ELEMENTS  WITH  ELEMENT  CONNECTIVITY  GIVEN  BY  MATRIX  CONN%%%%%%%%%% 

%%%%%%%%%%%% 

goblk=zeros(nunidof); 

7  5  gobIiii=zeros(nuindof); 

goblc=zeros(nuindof); 
for  i=l:nuniel 
v=comi(i,l); 
w=coiin(i,2); 

80  goblk(dof^node*  V- 1  :do0^ode*  v,dof_node*  v- 1  :dof_node*  v)= 

goblk(dof_node*  V- 1  :dof_node*  v,dof_node*  v-1  :dof_node*  v)+. . . 
ke(  1  :dof^node,  1  :dof_node); 

goblk(dofjiiode*v-l  :dof_node*  v,dof_node*w-l  :dof_node*w)=... 
goblk(dof_node*  V- 1  :dof_node*v,dof_node*w-l  :dof_node*w)+.. . 

85  ke(  1  :dof_node,dof_node+1 :2*dof_iiode); 

goblk(dof_node*  w- 1  :dof_uode*  w,dof_node*  v- 1  :doOiode*  v)= . , . 
goblk(dof^node*w-l  :doOiode*w,dof_node*v-l  :dof_node*v)+... 
ke(dof_uodei- 1 :2*dof_node,  1  :dof_node); 

90  goblk(dof_node*w- 1  ;dof_node*w,dof_node*\v-l  :dof_node*\v)=. . . 

goblk(dof_node*w-l  :dof_node*w,dof_node*w-l  :dof_node*w>f. 
ke(dof_node+l  :2*dof_node,doflnodef  1 :2*dof_node); 

goblm(dof_node*v-l:doOiode*v,dof_node*v-l:dof_node*v)=... 

95  gobliTi(dof_node*  V- 1  :dof jiode*  v,dof_node*v-l  :dof jiode*v>t-... 

me(  1  :dof_node,  1  :dof_node); 

gobLn[i(dof_node*  V- 1  :dof_node*v,dof_node*w-l  :dof_node*w)=.. . 
goblm(dof_node*  V- 1  ;dof_node*  v,dof_node*  w- 1  :dof_node*  w)-i-. . . 
nie(  1  :dof_node,dof_node+ 1 :2*dof_node); 

1 00  goblm(dof_node*w-l  :dof_node*w,dof_node*v-l:dof_node*v)=... 

goblin(dof_node*  w- 1  :dof_node*  w,dof_node*  v- 1  :dof_node*  v)+. 
nie(dof_node+l  :2*dof_node,  1  :dof_node); 
goblm(dof_node*w-l  :dof_node*w,dof_node*w-l  :dof_node*w)=... 
goblni(dof_node*w-l  :dofjiode*w,dof_uode*vv-l  :doOiode*w)+.. . 

1 05  me(dof_node+ 1 :2  *dof_node,dof_node+ 1 :2  *dof_node); 

end 

gobIc=sqrt(- 1  )*struc__danipmg.  *gobIk; 

goblkx=goblk; 

goblcx=goblc; 

110  goblnix=goblin; 


for  i=l  :iiuniel 
v=comi(i,l); 

115  \v=conn(i,2); 

goblkx(dof^node*  V- 1  :dof_node*  v,dof^node*  v-1  :dof^node*  v)= 
goblkx(dof_node*v-l  :dof^node*v,dof_iiode*v-l  :dof_node*v)+... 
Imnpspring(i).  *ke(  1  ;dofjiode,  1  :dof_node); 

1 20  goblkx(dof_iiode*  V- 1  :dof_node*  v,dof_node*w-l  :dof_node*\v)=. 
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goblkx(dof_iiode*  V- 1  :dof_ttode*  v,dof_ncMle*  w- 1  :dof_node*  w)f . . . 
limipspring(i).*ke(l  :dof_node,dof,node+l  :2*dof_node); 

goblkx(dof_node*w- 1  :dojf_node*  w^dof^node*  v- 1  :dof^node*  v)=. . . 
goblkx(dof^node*w-l  :dof_^node*w,do^iiode*v-l  :dof_node*v)4-... 
lmnpspring(i)*ke(dof_node+l:2*dof_node,l;doOiode); 

goblkx(dof_node*w-l  :dofjK)de*w,dof^node*w-l  :dof_node*w)=.,. 
goblkx(dof^node*w-l  :dof_node*w,dof_node*w-l  :do^ncxie*w)+... 
Imnpspring(i)  *ke(doOiode^-l  :2*dof_n(>de,dof.node+l  :2*doOiode); 

gobhiix(dof^iiode*v- 1  :dof_node*  v,dof_node*  V“  1  rdof^iiode*  v)=, . . 
goblnix(dof,node*v-.l  :dof_ncxie*v,dof_node*v-l  :dof_node*v>f ... 
Imiipmass(i).  ♦me{  1  :dof jiode,  1  :doOiode); 

goblinx(dofjttode*v-l  :dof^node*v,do^node*w-l  :dof^uode*w)=. .. 
goblinx(dof_node*v-l:dofjiode*v,dof_iiode*w-l:dofjiode^w>4'... 
lmiipniass(i).  *me(  1  :dofjiode,dof_node+ 1 :2*dofjiode); 

gobliiix(dof_node^w«'  1  :dof_node*w,dof_node*  v- 1  :dof_iiode*  v)=. . . 
goblinx(dof_iiode*w-l  :doOiode*w,dof_iiode*v-l  :dof_node*v>i-... 
lunipmass(i).  *me(dof_node+l  :2*dof_n^e,  1  :dof_node); 

gobliiix(dof_node*w-l  :dof_uode*  w,dof_iiode*  w- 1  :dof_node*w)=. . . 
goblmx(dof_node*w-l  :dof_node*w,dof_node*w-l  ;dof_node*w)4-. .. 
lmnpiiiass(i).  *me(dof_node+l  :2*dof_ncxie,dof_uodef  1 :2  *dof_node); 


goblcx(dof jiode*  V- 1  ;dof_node*  v,dof_node*  v- 1  ;dof_node*  v)=. . . 
goblcx(dof_node*v-l  :dof^node*v,do^node*v- 1  :dof_node*v)+. . . 
luinpdamp(i)*sqrt(-l  )*struc_damping,  *ke(  1  :dof_node,  1  :dof_node); 

goblcx(dof_node*v-l:dof_node*v,dof_node*w-l:dof_node*w)=... 
goblcx(dof_node*v- 1  :dof_node*  v,dof_node*w- 1  :dof_node*  w)+. . . 
luinpdainp(i)*sqrt(-l  )*struc_damping.  *ke(  1  :doOiode,doOiode^  1 :2  *doOiode); 

goblcx(dof_node*w-l  :dof_iiode*w,dof_node*v-l  :dof_node*v)=.. . 
goblcx(dof_node*w-l  :dof_node*w,dof_node*v-l  :dof^node*v)+... 
liimpdainp(i)*sqrt(-l  )*struc_dampmg.  *ke(dof_node+l  :2*dof_node,  1  :dofjiode); 

goblcx(dof_node*  w- 1  :dof_node*w,dof_node*  w- 1  :dof_node*w)=. . . 
goblcx(dof_node*w-l  :dof_node*vv,dof_node*w-l  :dof_ncxle*w>f-... 
lumpdainp(i)*sqrt(-l  )*struc_daiiiping.  *ke(dof_node+l  :2*dof_node,dof_node+l  :2*dof_node); 
end 


iiuiudof=nuindof-doftokill; 

[gk,gin,gc,k_ai\al,m_anal,c_anaI]=fIxbcs(goblk,goblm,gobIc,bc); 

[gkx,gmx,gcx,k_exp,m_exp,c_exp]=fixbcs(goblkx,goblmx,goblcx,bc); 
save  beaindata 
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FSTATIC.M 


function  [kstat,mstat]=fstatic(k4n,oset,aset) 

aset__size=lengtli(aset); 

kaa=k(aset,aset); 

kao=k(aset,oset); 

5  koo=k(oset,oset); 

koa=kao’; 
clear  k; 

k=[kaa,kao,koa,koo]; 

1 0  maa=in(aset,aset); 

mao=m(aset,oset); 
moo=m(oset,oset); 
moa=mao'; 
clear  m; 

1 5  m=[inaa,mao;moa,moo]; 


t_static=-koo\koa; 
T_static=[eye(aset_si2e);  t_static]; 

20 

kstat=T_static’*k*T_static; 

.  nistat=T_static’*m*T_static; 

end 

25 
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FIRS  TAM.M 


% 

function  [kirs,niirs]=firs_taiii(k,ni,oset,aset) 

% 

%  this  function  returns  tlie  IRS  reduced  stiffness 
%  and  mass  matrices,  given  tiie  unreduced  couterparts. 

%  Care  must  be  taken  tliat  tlie  aset  and  oset  vectors  correspond 
%  witli  tlie  existing  arrangement  of  k  and  m. 

%  k  and  m  are  UNPARTITIONED  matrices 
% 

aset_size=lengtli(aset); 

% 

kaa=k(aset,aset); 
kao=k(aset,oset); 
koo=k(oset,oset); 
koa=kao'; 
clear  k; 

k=[koo,koa;kao4caa]; 

% 

maa=m(aset,aset); 
mao=m(aset,oset); 
moo=m(oset,oset); 
moa=mao’; 
clear  m; 

m=[moo,moa;mao,maa]; 

% 

t_static=-koo\koa; 

T_static  =  [t__static;  eye(aset_size)]; 

% 

kstat=T_static’*k*T_static; 

mstat“T_static’*m*T_static; 

% 

^i^s=t_static+inv(koo)*(moa+moo*t_static)*mv(mstat)*kstat; 

T_irs=[tirs;eye(aset_size)]; 

% 

kirs=T_irs'*k*T_irs; 

rairs=T_irs’*m*T_irs; 

% 

%  end  function  firs_tam 
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FREQMODE.M 


%%%%%%%%%%%%%%%%%%%%%0/o%%0/o0/o0/o%0/o0/^,0/^,0^0/^0/o0/^0/^0/^,0/^0/^0/^0/^0^0/o0^^^^ 

%%%%%%%%%%% 

%tliis  function  returns  a  vector  U  containing  modal%%%%%%%%%%%%%%%%%%%%% 

%frequencies  (rad/sec)^2  in  ascending  order  along%%%%%%%%%%%%%%%%%%%%%% 

5  %witli  associated  mode  sliapes%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%o/o%%%o/oO/oO/oO/oO/^o/^ 

%%%%%%%%%%%%%%%%%%%%%%%%%o/o%%%%%%%o/o%%%o/o%%%%%%%%%o/oO/oO/^,o^o/^o/^,o^o/^o^o/^o/^^^^^^ 

%%%%%%%%%%% 


1 0  function  [ul , lambda, index]=ffeqmode(k,m) 

[u,lambdat]=eig(in\k); 
[lambda,index]=sort(diag(iambdat)); 
ll==zeros(length(k),length(k)); 
for  i=l:iengtli(k); 

15  ll(i,i)=lambda(i); 
end 

Iambda=diag(li); 
forj  =  l:length(k) 
ul(:J)  =  u(:,index(j)); 

20  end 
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NDX3D.M 


5 


10 


15 


bo 


function  [r,c]  =  ndx3d(siz,ij,k) 

%NDX3D  Index  into  3-D  matrix  packed  in  a  2-D  matrix. 

%  ELEM  =  NDX3D([M  N  P],I,J,K)  returns  tire  element  position  ELEM 
%  of  tlie  (i  j,k)  elements  of  a  M-by-N-by-P  matrix  wltich  is 

%  stored  in  a  (M*N)-by-P  matrix.  For  example,  tlie  tliree  M-by-N 
%  matrices  A1  ^2^3  are  packed  into  a  2-D  matrix  using 
%  A=[A1(:)A2(:)A3(:)]; 

%  If  lengtli(I)  is  m,  LENGTH(J)  is  n,  and  LENGTH(K)  is  p, 

%  then  ELEM  will  be  an  m-by-n-by-p  matrix 
% 


%  t^>9l  ~  NDX3D([M  N  P],I,J,K)  returns  tlie  row  and  colunm 

%  position  of  the  (i  j,k)  element  as  stored  in  tlie  normal 

%  2-D  matrix  of  size  (M*N>by-P. 

% 

%  To  specify  all  the  elements  along  one  dimension  use 
%  For  instance,  NDX3D([3  5  4],2:3,':',3;4)  returns  tlie 
%  elements  for  the  2-by-5-by-2  matrix. 

% 


%  See  also  ELEM3D,  MESHGRID,  SLICE. 


%  Clay  M.  Tliompson  1 1-3-92 
%  Copyriglit  (c)  1992  by  Die  MatliWorks,  Inc. 

%  SRevision:  1.7$  SDate:  1993/09/03  14:36:52$ 


r5  if  isstr(i),  i  =  1  :siz(l );  end 

ifisstr(j)J  =  l:siz(2);end 
if  isstr(k),  k  =  1  :siz(3);  end 


if  isempty(i)  |  isenipty(j)  |  isempty(k),  r  =  [];  c  =  [];  return,  end 

10 

if  any(  (i<0)  |  (i>si2(l)) ),  errorChidex  I  out  of  range.');  end 
if  any(  O'<0)  |  (j>siz(2)) ),  errorChidex  J  out  of  range.');  end 
if  any(  ^<0)  |  (k>si2(3)) ),  eiTor(1ndex  K  out  of  range.');  end 

5  ifnargout==2 

Gj,ii]  =  meshgrid(j,i); 
r  =  ii(:)  +  (jj(:).l)*siz(l); 
c  =  k(:); 


0  else 

[ij,iijkk]  =  nieshgrid(j,i,k); 
r  =  ii  +  (jj-1  )’»siz(  1 )  +  (kk-1  )*prod(siz(  1 :2)); 

end 
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INTSUB.M 


fimction  [K,M,C]-intsub(Z__I,  Z_R,  omega,  omega  1,  W,  set^size) 
waitbar_haiidIe=waitbar(0, ’Computing  Integrals’); 
j=sqrt(-l); 


5  mteg=omega.'’^2.  *abs(W); 

rint=real(integ); 
iint=imag(mteg); 

aa=mytrapz(omegal  ^t)fmytrapz(omegaI  ,iiut)*j; 
mteg=(L/omega.^2).*abs(W); 

1 0  rint=real(integ); 

imt=imag(mteg); 

bb=mytrapz(omegal  ,rint>{-mytrapz(omegal  ,iiut)*j; 

integ=abs(W); 

rmt=real(integ); 

15  iint=imag(mteg); 

cc=mytrapz(omegal,rmt>i-mytrapz(omega  1  ,iiiit)*j; 
for  i=l  :set_size 
waitbar(i/set_size); 
fork=i:set_size 

20  M(i,k)=(  1  /(aa*bb-cc'^2  ))*(bb*mytrapz(omega  1 . 

omega.  ♦ZJ(ndx3d([set_size  set_size  lengtlXomega)],... 
i,k,  1  :leiigtli(omega))).  ’*‘abs(  W))-cc*mytrapz(omega  1 . 

1, /omega  *ZJ(ndx3d([set_size  set_size  lengtli(omega)]... 
,i,k,l:length(omega))).*abs(W))) ; 

25  M(k,i)=M(i,k); 

C(i,kH  l/cc)*mytrapz(omega  1  ^_R(. . . 
ndx3d([set_size  set_size  lengtlTi(omega)],... 
i,k,  1  :lengtli(omega))).  *abs(W)); 

C(k,i)=C(i,k); 

30  K(i,k)=(  l/(aa*btKX^2))*(cc*mytrapz(omega  1 . 

omega.  *ZJ(ndx3d([set_size  set_size  length(omega)]... 
,i,k,  1  :length(omega))).*abs(W))-aa*mytrapz(omega  1 ,. . . 
l./omega.*ZJ(ndx3d([set_size  set_size  length(omega)]... 
,i,k,l:lengtli(omega))).*abs(W))) ; 

35  K(k,i)=K(i4c); 

end 
end 

close(waitbar_liandle) 
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MYTRAPZ.M 


function  z  =  trapz(x,y) 

%TRAPZ  Trapezoidal  numerical  integration. 

%  Z  =  TRAPZ(X,Y)  computes  tlie  integral  of  Y  witli  respect  to  X  using 
%  trapezoidal  integration.  X  and  Y  must  be  vectors  of  tlie  same  lengtli, 
%  or  X  must  be  a  coluima  vector  and  Y  a  matrix  witli  as  many  rows  as  X. 
%  TRAPZ  computes  tlie  integral  of  each  column  of  Y  separately. 

%  Tlie  resulting  Z  is  a  scalar  or  a  row  vector. 

% 

%  Z  =  TRAPZ(Y)  computes  tlie  trapezoidal  integral  of  Y  assuming  unit 
%  spacing  between  tlie  data  points.  To  compute  tlie  integral  for 
%  spacing  different  from  one,  multiply  Z  by  tlie  spacing  increment. 

% 

%  See  also  SUM,  CUMSUM. 

%  Clay  M.  Tliompson,  10/16/90;  Cleve  Moler,  1/19/92. 

%  Copyriglit  (c)  1984-94  by  Tlie  Matli Works,  Lie. 

%  Make  sure  x  and  y  are  column  vectors,  or  y  is  a  matrix. 

%  Trapezoid  sum  computed  witli  vector-matrix  multiply. 

[ni,n]=size(x)  ; 
z=0  ; 

for  mdex=l:m 

y’y=y(  l+(ui(iex-l  )*n:n+(mdex-l  )*n); 

yy=yyCX 

xx=x(iiidex,:); 

XX  =  xx(:); 

z  =z+difr(xx)’  *(yy(l:n-l)  +  yy(2:n))/2; 
end 
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